9 ' 




OF TECHNOLOGY 

INVESTIGATION OF THE 
EPOCH STATE FiLTER 

by 

Joan Annette Edwards 

February 1972 



PREPARED AT 

CHARLES STARK DRAPER LABORATORY 

CAMBRIDGE, MASSACHUSETTS, 02139 




Reproduced by 

NATIONAL TECHNICAL 
INFORMATION SERVICE 

US Department of Commerce 
Springfield, VA. 22151 







INVESTIGATION OF THE EPOCH STATE FILTER 


by 

JOAN ANNETTE EDWARDS 
B. S. A. E. , PURDUE UNIVERSITY 
(1970) 


SUBMITTED IN PARTIAL FULFILLMENT 
OF THE REQUIREMENTS FOR THE 
DEGREE OF MASTER OF SCIENCE 
at the 

MASSACHUSETTS INSTITUTE OF TECHNOLOGY 


February, 1972 


Signature of Author 

Certified by 
Accepted by 



_Q C AvUJTLArjoD 


Department of Aeronautics and 
Astronautics, January 21, 1972 

'll 

Thesis Supervisor 



Chairman, Departmental j 
Graduate Committee 


I 



ACKNOWLEDGMENT 


I would like to express my appreciation to Dr. Richard H. 

Battin, my thesis supervisor, for suggesting that I develop as a 
thesis subject. The Epoch State Filter, which he initially conceived 
and for reviewing my final paper. I especially owe a debt of gratitude 
to Dr. Steven R. Croopnick for his continued assistance and encourage- 
ment through my thesis work. 

I would further like to thank the Control and Flight Dynamics 
Division of the Charles Stark Draper Laboratory for providing this 
research opportunity and allowing me to use its research facilities. 

A special note of appreciation is extended to the National 
Aeronautics and Space Administration and to ZONTA International 
for sponsoring my graduate study at M. I. T. 

This report was prepared under DSR Project 55-41200, 
sponsored by the Manned Spacecraft Center of the National Aero- 
nautics and Space Administration through contract NAS 9-10386. 

The publication of this report does not constitute approval 
of the National Aeronautics and Space Administration or the Charles 
Stark Draper Laboratory of the findings or the conclusions contained 
herein. It is published only for the exchange and stimulation of ideas. 


ii 




INVESTIGATION OF THE EPOCH STATE FILTER 


by 

Joan Annette Edwards 

Submitted to the Department of Aeronautics and Astro- 
nautics on January 2l, 1972 in partial fulfillment of the require- 
ments for the degree of Master of Science. 

ABSTRACT 


A new navigation filtering technique has been formulated using 
as state variables the initial or epoch position and velocity of the 
spacecraft. The estimate of this initial state is then improved 
by filtering new measurements. The current state may be obtained 
by a conic extrapolation of the epoch state. Results of a digital 
computer simulation of the epoch state filter show that this formu- 
lation of the navigational problem results in less computer run 
time and less computer storage space than conventional techniques. 
The errors produced by this technique have been demonstrated to be 
comparable to those obtained by conventional maximum- likelihood 
filtering. 
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LIST OF SYMBOLS 


General Notation 

An underlined symbol indicates a vector. 

A prime to the upper right of a symbol indicates the 
quantity is that extrapolated from the previous 
measurement time. 

A caret over a symbol indicates that the quantity is an 
estimate. 

A bar over a symbol or group of symbols indicates the 
expected value of what is beneath. 

A "o" subscript on a symbol denotes an epoch quantity. 

A "k" subscript on a symbol denotes that quantity at the 
time of the measurement. 
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disturbing acceleration 
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apriori variance of measurement error 
geometry vector 
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T- 1 T 
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variable defined in t^ equation 
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vii 




6 

e 


F,(t) 


°c<‘> 

°ad<‘> 


variation in true anomaly difference 0 
position deviation from osculating orbit 

error vector 

error incurred using ESF 
error in position estimate 
extrapolated covariance matrix 
covariance matrix of estimation errors 
extrapolated epoch covariance matrix 
epoch covariance matrix 
conic epoch covariance matrix 
epoch state filter 

variation in generalized anomaly x 
true anomaly 

scalar quantity in Y equation 
scalar quantity in Y equation 
matrix in ^ equation 
conic F(t) 

gravitational force per unit mass 
scalar quantity in Y equation 
scalar quantity in Y equation 
gradient of gravity vector 
conic G(t) 

G(t) - G^(t) 

angular momentum 

unit vector in radial direction 

unit vector normal to i^ in orbital plane 

viii 




1 

— z 


I 

J. 

R 


6q 


6q' 


E 


r 

— o 


r 

-osc 

R 

R* 

t 

o 

U 

u 


n 


V 

V 

V 




V 

-osc 


i X 

identity matrix 
second zonal harmonic 
error in velocity estimate 
orbital parameter 
Encke integration variable 
measured quantity 

extrapolated estimate of measured quantity 

position vector 

position 

equatorial radius of earth 
epoch position vector 
epoch position 
osculating position vector 
ar/avo 

ar^/av 

epoch time 

gravitational parameter 

n*''^ transcendental function 

velocity deviation from osculating orbit 

velocity vector 

velocity 

epoch velocity vector 
epoch velocity 
osculating velocity 
av/av 


IX 



V* 

i£ 

— o 
X 


6x 

6x 
— o 

Y 

a 

0 

60 

A0 

cos 0 

$ 


c 




weighting vector 
epoch weighting vector 
generalized anomaly 

state vector 

estimated state deviation vector 

estimated epoch state deviation vector 

epoch to current state transformation matrix 

variational orbital element 

variational epoch orbital element 

true anomaly difference 

deviation in true anomaly difference 

epoch 0 deviation 

total change in 0 

i • i 
— r — z 

state treinsition matrix 
conic state transition matrix 



X 



CHAPTER I 


INTRODUCTION 


One of the processes of coasting flight navigation involves 
improving the estimate of the spacecraft's position and velocity 
vectors. In the Apollo navigation system this is accomplished 
with a recursive formulation of the maximum- likelihood es- 
timator in which state or process noise is neglected. Measurement 
data, regarded as scalar information, is incorporated as it is 
obtained to update the estimate of the state. This data handling 
technique is termed recursive processing as opposed to batch 
processing wherein measurement data is incorporated all at once. 


The process of determining the estimate of the state vector 
as stored in the on-board digital computer involves integrating 
the equations of motion which govern the spacecraft. During 
both coasting flight and under the influence of a vector disturbing 
acceleration, a^, these equations are: 


d^r(t) 


d t 


2 


^r(t) 

r^(t) 


a^(r(t)) 


dr(t) 

v(t) = = 

dt 


where r(t) and v(t) are the current position and velocity vectors 
of the vehicle with respect to the primary body and ju is the 
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gravitational constant of the primary body. Integration of these 
equations involves selecting an appropriate set of state variables. 
For Apollo, current position and velocity are used to define the 
state. Intuitively, this is a proper choice of state variables since 
cxirrent position and velocity are ultimately the quantities which 
are estimated by the navigation filter. Another related set of 
variational parameters is the spacecraft's position and velocity 
at some initial time or epoch. 

The formulation of the space navigation problem using the 
epoch state as state variables requires a variation of parameters 
solution which is discussed in this thesis. It is shown that this 
latter formulation has certain computational advantages over 
the conventional one, namely, less computer run time and computer 
storage space. The disadvantage, a slight decrease in accuracy, is 
introduced because of a simplifying assumption used to integrate 
the time derivatives of the epoch state error covariance matrix. 
However, in many cases this error is small so that the epoch 
formulation of the navigational problem may replace a formulation 
such as the one used in Apollo. Various examples are given as 
cases where the implementation error is negligible. Also, 
statistical equations were developed to predict the filtering 
approximations of the epoch state filter. 

In Chapter II the conventional formulation of the naviga- 
tional problem is discussed in detail. Its design philosophy is 
explained as well as how error is propogated and measurement 
date is recursively incorporated to improve the estimate of the 
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state vector. Also discussed in this chapter is the design philo- 
sophy and choice of state variables for the epoch formulation 
of the navigational problem. The equations for the epoch filter 
used to estimate the state vector are derived in Chapter III. 
Explanation is given for the assumption made to simplify this 
formulation. A comparison between the basic equations of the 
conventional and epoch filter formulations is also given. Finally, 
statistical equations for the error incurred using the epoch 
rather than the conventional state filter are derived. Computer 
simulation results of the epoch formulation of the navigational 
problem are presented and discussed in Chapter IV. Conclusions 
regarding the epoch state filter, its advantages and areas of 
application are explained in Chapter V, 
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CHAPTER II 


NAVIGATION FILTER FORMULATIONS 

2. 1 Design Philosophy for the Apollo Navigation Filter 

Position and velocity as maintained in the Apollo Guidance 
Computer (AGC) are estimates of the true state of the spacecraft. 

These estimates are propagated from measurement to measurement 
by integrating the equations of motion of the vehicle with respect to 
time. Integration of the spacecraft's motion involves the selection 
of an appropriate set of state variables. Current position and velocity 
of the spacecraft are used in the Apollo navigation system. Intuitively, 
this is a proper choice of state variables since these are the quantities 
to be estimated. 

The vector equations governing the motion of the spacecraft 
during coasting flight are: 

d^r(t) 

— + ^ r (t) = a. (r(t)) (2.1.1) 

dt2 r3(t) - “ 

dr(t) 

= V (t) (2.1.2) 

dt ~ 

where r (t) and v (t) are the vector position and velocity of the vehicle 
in non- rotating rectangular coordinates with respect to the primary 
body and /x is the gravitational parameter of this body. The 
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quantity is the vector disturbing acceleration which prevents 

the motion of the spacecraft from being precisely a conic orbit. The 

disturbing acceleration is a function solely of the position vector. 

For earth orbit, only gravitational perturbations due to the non- 

spherical gravity field of the earth need be considered in a The 

— d 

equation for the disturbing acceleration used in this study is given by: 

r 2 

id "" - ~ ^2 (~) [ - 5 cos0)i + 2 cos 0 i ] (2. 1. 3) 

2 ^ ^ ^ " 

where 0 is the angle between _i^ the unit vector in the r direction, and 
_i^, the unit vector in the direction of the spin axis; r^, is the equatorial 
radius of the earth and J2 is the coefficient of the second harmonic of 
the earth's potential function. 

When is small compared with the central field of the primary 
body, direct integration of Eqs. 2. 1. 1 and 2. 1, 2 in rectangular coordi- 
nates is inefficient. An alternate procedure suggested by Encke^, is 
used to perform this integration. For the Encke method of integration, 
the actual position and velocity of the spacecraft, defined by the current 
values of r(t) and Wt), are viewed as deviations from a conic or oscu- 
lating orbit 

= lose ^ (2. 1.4) 

" Zosc (2.1.5) 

In practice, the osculating orbit and the deviations from this orbit 
onboard a spacecraft are only estimates of their true values and are 
represented with a superscript "" ". Hence, the current position and 
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velocity estimates are given by 


r(t) = r^^^(t) 6 


v(t) - V (t) + ^(t) 
— — ” O S C ““ 


The osculating orbit at any particular time is determined from 
ideal two-body motion by solving Kepler's equation for 9. This is 
accomplished by using the following equations for two-body motion: 


r = F r + G V 
— osc — o — 


where 


V = F.r + G, V 
— osc t — o t — o 


1 - — (1 - cos 0) 
p 


r — sin 9 
h 




— (1 - cos 9) - -i— sin 0 
r p r h 


G, = 1 (1 - cos 9) 

t p 


1 + (— - 1 ) cos 9 — ^ — 2 . s in 9 
^o ^o 


The parameters in 2, 1. 14 are determined according to the 
following equations: 
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r 

— o 


V 

— o 


(2, 1. 15) 


cr 


o 




Oi 




(2. 1. 16) 


P " 



(2. 1. 17) 


The deviation vector ^(t) and y(t) are obtained by integrating the 
following differential equations: 


dy(t) 

dt 


d6(t) 

= V (t) 

dt ~ 


r-" (t) 

osc 


-[ f(q) r(t) +_6(t)] +^(r(t)) 


(2. 1. 18) 
(2. 1. 19) 


subject to initial conditions = jf(tj^) = 0 where 


and 


„ _ i • (A - 2x) 

^ O 


f (q) = 


q(3 + 3q + q^) 

1 + (1 +q)^/^ 


( 2 . 1 . 20 ) 


( 2 . 1 . 21 ) 


A recommended numerical integration technique, Nystrom’s 

Method, exploits the fact that a^ is a function only of r, the vector to 
be integrated. 

For Encke's Method to be efficient, the first term on the right 
hand side of Eq, 2. 1,19 must remain small, i, e, , of the same order or 
less as the disturbing acceleration. To insure the efficiency, a new 
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osculating orbit is periodically defined from which ^ and^ are cal- 
culated. When this rectification is done, the new osculating orbit 
is defined by the current values of r (t) and v(t) and the initial condi- 
tions for Equations 2. 1. 18 and 2.1 . 19 are again set equal to zero at 
the current time. 


2.2. Error Covariance Matrix 


The position and velocity vectors which are stored in the Apollo 
Guidance Computer (AGC) are estimates of their true values. Since 
these estimates will be in error, it is necessary as part of the 
maximum- likelihood filtering technique to maintain the statistics 
associated with these errors. 

If € (t) is the three dimensional error in the position estimate 
and T)(t) is the three dimensional error in the velocity vector, then 
the error in the estimate of the state vector is given by 


e (t) 


_f (t) 
THt) 


( 2 . 2 . 1 ) 


When unbiased measurement data is processed to determine 
the maximum - likelihood estimate, the error in the estimate has a 
zero mean, i, e. e = 0 so that the 6x6 covariance matrix of estima- 
tion errors is defined by 


T 

E = e e ^ = 


T 

c e e T1 


T T 

n C " 7] Tl-" 


( 2 . 2 . 2 ) 
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and is also stored in the AGC. 


A useful measure of the error in the position estimate is given 
by the rms position error. This error is determined by the square 
root of the trace of the upper left hand 3x3 partition of the covariance 
matrix and is given by 


rms position error = 



1/2 


( 2 . 2 . 3 ) 


Similarly, the rms velocity error, a good measure of the error in the 
velocity estimate, is determined from the lower right hand corner 
of the covariance matrix according to 


rms velocity error = 



( 2 . 2 . 4 ) 


With the recursive formulation of the Kalman estimator, mea- 
surement date is processed as it is obtained. The covariance matrix 
is maintained in the AGC in the intervals between which measurements 
are taken and is updated as is the current estimate of the state vector 
when the measurement data is incorporated by this linear estimator. 

The Kalman filter operates as follows. First, the old estimate 
is extrapolated to the current time, yielding the best estimate prior 
to the incorporation of measurement data. For coasting flight 
E (tj^ j) is extrapolated to the current time, t^^, by 


*<‘k- tk-/ 


( 2 . 2 . 5 ) 
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The prime ' to the upper right of E(tj^) indicates the covariance matrix 
of estimation errors at t^^ is that based on previous k-1 measurements 
and $ (t, , t, . ) is the 6x6 state transition matrix by which the 

cC K. " i 

state and certain statistical quantities are extrapolated in time from 

t to t, . The transition matrix satisfies the first order matrix 
k-1 k 

differential equation 

^<^k *k-i) = F(t) tj^_ j) (2.2.6) 

subject to the initial condition 


* <‘o' *0> = ' 


where I is the 6x6 identity matrix. 


F (t) 


0 I 


G(t) 0 


where 


G(t) = \\— II 

9r 


(2. 2. 7) 


The 3x3 matrix G(t) is the gradient of the gravitational field ^ with 
respect to the components of the position vector r. For orbital 
navigation about a primary body, G(t) is given by 
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G(t) 


( 2 . 2 . 8 ) 


= — [ 3 r (t) r (t)'^ - (t) l] 

r5(t) “ - 

An alternate method of extrapolating the covariance matrix 
rather than by first determining j) and then substituting 

this matrix into Equation 2.2.5 is to integrate the first order differen- 
tial equation for E'(t ) 

K. 

E'(tj^) = F(tj^) E'(t^) + E'(t^) F(tj^)'^ (2.2.9) 

This is obtained by differentiating Equation 2. 2.5 with respect to time 
and substituting Equation 2. 2.6 in the resulting equation, the deriva- 
tion of which is given in Appendix A. 

Once the extrapolated covariance matrix is obtained, the 
measurement data is incorporated according to optimal estimation 
theory. As a result of the measurement incorporation, the statistics 
of the error covariance matrix are changed. A weighting vector co 
is determined which minimizes the mean squared error in the esti- 
mate. According to maximum-likelihood theory, the weighting 
vector is given by 

CO = J- E'b (2. 2. 10) 

— a — 

where b is a 6 dimensional geometry vector associated with the 
measurement and 

a = b”^ E'b + (2. 2. 11) 
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where a is apriori variance of the measurement error. In terms 
of E'(t, ) as well as Equations 2. 2,10 and 2. 2, 11, the new value for 
the covariance matrix of estimation errors is determined according 
to 


E (t^^) = E’(tj^) - w e' (t^) 


( 2 . 2 . 12 ) 


or 


E(tj^) = (I - w b'^) e' (tj^) 


(2. 2. 13) 


2.3 Measurement Incorporation 


For flight paths which are close to a nominal one, 6 x may 
be expressed as a linearized deviation about the nominal state and 
is denoted by 


6x 


6 r 
6v 


(2.3. 1) 


The estimate of the state vector is obtained by the operation 
of the optimum linear estimator on the state deviation vector. First, 
the previous state deviation estimate is extrapolated to the current 
time yielding its best value prior to the incorporation of new infor- 
mation. This is expressed by the following relationship 


6x'(tj^) = j) ^(tj^_j) 


(2. 3. 2) 
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V 

where $ (tj^, j) is the state transition matrix. The best estimate 
of the measured quantity 6q', is computed according to 


6q'(tk) = 6x' (t^) 


(2. 3. 3) 


The difference between the actual measurement data 6q and 
the filter's prediction of what this value should be 6q' is weighed 
statistically against 6x'. This is accomplished through the use 
of a statistical weighting vector to defined by 


E'b 

W = (2.3.4) 

a 


where 


a 


T ? 

b E ' b + 01 


(2. 3. 5) 


Making use of Equations 2. 3. 2 through 2. 3. 5, the updated state 
estimate at measurement time t^ is obtained from 


6x(tk) = 6x/ (t^) + w(6q (tk) - 6q' (tk)) 


(2. 3. 6) 
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Equation 2. 3. 6 is simplified by adding the estimate of the state 


deviation from the nominal path to the current state; 


2nom<‘k> ' i<‘k> *2.3,7) 

SO that a new nominal path is defined at every measurement time t^ . 
By adding 6x to the current estimate of the state, the spacecraft 
is assumed to be on the nominal trajectory which is redefined at 
each measurement time. The extrapolated deviation of the vehicle 
from the nominal path at time tj^ is then zero since the spacecraft 
was on the nominal trajectory at time tj^ This is illustrated in 
Figure 2.1. 

The equation for determining the state deviation estimate 

at t, then reduces to 
k 


fix (tj^) = w fiq (tj^) (2.3.8) 

A I 

This is seen by substituting fix (t ) = 0 in Equation 2. 3. 2 and noting 

— K — 

that 6 q' which is a function of 6x' (t ) is also zero. 

K 

2. 4 Design Philosophy of the Epoch State Filter 

Current position and velocity of the spacecraft are the typical 
vector quantities which are estimated by a spacecraft computer doing 
recursive processing such as the AGC. These quantities vary con- 
tinuously along the path. An estimate of the current state is obtained 
by integrating the equations governing the motion of the vehicle. In 
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Figure 2. 1 Effect of adding the expected deviation of the 
state S^to the expected state i at' time t|< 
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the Apollo formulation of the navigational problem, current position 
and velocity along the path of the spacecraft are the quantities 
used in the integration of these equations. This is a convenient 
and intuitively correct choice of state variables since current position 
and velocity are the quantities to be estimated. Estimation of 
the state vector is then made by incorporating measurement data 
using Kalman gains in the space navigation filter. This current 
measurement information is used to update the current estimate 
of the position and velocity vectors. 

The formulation of the navigational problem developed in 
this paper employs a related set of state variables, namely, the 
position and velocity of the spacecraft at some initial or epoch 
time which are adjusted in a manner such that a simple conic 
trajectory connects this point with the current position. This 
adjustment is made so that the extrapolated velocity also matches 
the current velocity at the current time. A scalar variational 
parameter is also used along with the above two vector quantities, 
as a means of extrapolating quantities from the epoch to the current 
state. This is the true anomoly difference, 0, the central angle 
between the current and the epoch position vectors and is considered 
the independent variable. 

Because the epoch state is not accurately known initially, 
measurement data is incorporated using a Kalman filter to estimate 
these initial conditions. This is analogous to what is done in Apollo 
however, for that formulation of the navigation problem, measurement 
data is used to improve the current estimate of the state. When 
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sufficient new information is incorporated, the epoch state vector 
is brought up to the current time by solving Kepler's equation. 

The main difference between the conventional and epoch 
formulation of the navigational problem is illustrated in Figure 2.2. 

In the conventional state filter, current position and velocity are 
used as state variables whereas for the epoch state filter, epoch 
position and velocity are used as state variables. These epoch 
quantities are integrated between measurements and then updated 
at the current time using current measurement data. 

The epoch formulation of the navigation problem developed 

2 

in this thesis makes use of the variable epoch form of the variational 
equations. This means that the epoch time, t^, is forced to vary 
in the intervals between measurements. The variable epoch form 
of these equations was used because of their relative simplicity 
when compared with the fixed epoch form, however, the navigational 
problem could just as well have been formulated using the fixed epoch* 
form of the variational equations for which t^ =0. This feature is 
explained in Chapter III. 

Rather than using the true anomaly difference 0 as the indepen- 
dent variable, the generalized anomaly, x, may be treated as the 
scalar variational parameter. Solving a differential equation for x 
eliminates the necessity of solving Kepler's equation for this same 
variable in the same fashion as integrating the differential equation 
for 0 eliminates obtaining Kepler's solution. Still another differential 
equation in terms of the epoch time, t^, can be integrated instead 
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Path of spacecraft 





Measurement 

Interval 


), E (t ) 
-0 0 0 0 


I 

I 






X(tN> 


Epoch formulation of the Navigational Problem 
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of the two equations indicated above. However, solving this 
equation for t^ necessitates solving Kepler's equation for x. 

Position and velocity were chosen as state variables for both 
the Apollo and the epoch formulations of the navigation problem for 
convenience, because this state is the quantity that is estimated by 
the filter. However, other variational parameters may be used 
to formulate estimates of the current state vector, such as the 
orbital elements. Equations for this formulation are developed 
in Reference 2, 
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CHAPTER III 


THE EPOCH STATE FILTER 


3. 1 Derivation of the Epoch State Filter 

Current position and velocity are the state variables used in 
the Apollo navigation filter and in the integration of the equations 
of motion of the spacecraft. The current state is estimated by 
measurement data which is incorporated to update the current state 
estimate directly. 

For the Epoch State Filter (ESF), position and velocity of the 
vehicle at the epoch time are the state variables which are used 
in the integration of the equations of motion. This filter estimates 
the current state indirectly by first estimating the epoch state. 
Current measurement data is processed to update the epoch state 
estimate. The improved current state estimate is then obtained 
from the epoch state estimate by conic extrapolation. 

The state equations of motion used for the ESF formulation 
are derived as variations of the epoch state, r^ and v^. An 
additional differential equation in terms of the independent scalar 
variable 0 is also integrated. Theta, the true anomaly difference, 
is the angle between the epoch and current position vector as 
illustrated in Figure 3.1. 
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These differential equations for r v and 6 (derived in 

—o' —o' 

Section 3. 3) are integrated to propagate the estimate of the epoch 

state before the new information is introduced. Measurement 

data is then incorporated to estimate the epoch state deviation, 

6 X (t). from its nominal value at the current time. This is related 
— o 

to the estimate of the current state deviation 6 x (t) by 

6 X = $ ^ 6 X (3. 1. 1) 

— o — 

in terms of the state transition matrix $ (t, t^). The current state 
estimate, 6x, determined by Equation 2. 3. 6 is rewritten here as 

6 X = 6 X + (jO (6q - b'^ 6 x' ) (3.1.2) 

replacing 6 q' by its equivalent b 6x where 

w = ^ (3.1.3) 

a 

The current covariance matrix of estimation errors, e', before 
the measurement is incorporated is related to the epoch covariance 
matrix by 


I I T 

E = E $ 
o 


(3. 1.4) 


Substituting Equations 3. 1. 2 through 3. 1. 4 into Equation 3. 1. 1 
and bracketing significant terms yields 
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6x = [ $ ^ 6x '] + 4>~ ^ 4>E ' /a) (6q - b'^ 6x ') (3. 1. 5) 

The product ^ $ is the identity matrix so that the first term in 

brackets is just 6 x^, i. e. 

6x' = 6x' (3. 1. 6) 

— — o 

and b is related to b according to 

b = b (3.1.7) 

— o — 

so the second term in brackets is b^. When Equations 3. 1. 6 and 
3. 1. 7 are substituted into Equation 3. 1. 5 the resulting equation is 

6Xq = 6x^' + [ E^' k^/a] ( Sq' - b”^ 6 X ' ) (3.1.8) 


The term in brackets can be defined as the epoch filter gain 
where 



(3. 1. 9) 


T 

Substituting Equation 3.1.9 into 3. 1. 8, post multiplying b in 
E:iuation 3. 1. 8 by the identity matrix and bracketing significant 

terms yields 


6x^ = 6 Xq' + Wq (6q- [b'^«ii] ^ 6x'] ) 


(3.1. 10) 


23 



But 


b'^ ^ = b (3.1. 11) 

— — o 

and 

= 6xJ (3.1.12) 

When Equations 3.1. n and 3. 1. 12 are substituted into Equation 3. 1. 10 
the resulting equation for the ESF is 

6x = 6 X ' + 00 (6q - b ') (3.1. 13) 

— o — O — O — O — O 

As is done in the Apollo navigation system, 6x^' is added to the 
estimate of the total state vector so that the nominal path is redefined 
by the measurement. Therefore, for the nominal path defined 
by the (k- 1)^*^ measurement, the extrapolated estimate of the 
deviation from the nominal trajectory is zero, i. e, 6x^' (tj^) = 0 
as illustrated by Figure 3,2. Substituting 6x^ = ^ into Equation 
3.1. 13 yields 



^o 


(3.1. 14) 


which is the ESF equation for incorporating scalar measurement 
data to update the estimate of the epoch state. 
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S) 



Figure 3. 2 Redefinition of the Nominal Path 
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3. 2 Derivation of the Error Covariance Matrix Used in the ESF 


The equation for the epoch error covariance matrix E^ 
associated with the state deviation 6x^ is obtained from the covariance 
matrix equation developed in Section 2, 2 


= "=k'- 


(3.2. 1) 


where the subscript k indicates the quantity at time tj^. Noting the 
relations between the epoch” and current parameters 


I • * T 

= #,..^E. * 


^k’O ^o 0 


(3. 2. 2) 


\ ^ %’O^o^k, 0 


(3. 2. 2) 


and substituting these into Equation 3. 2. 1 yields 




rp b b rp 

^ ^ - 4>E 

o o 


(3. 2. 4) 


where co in Equation 3. 2. 1 has been replaced by 



(3. 2. 5) 


and the subscripts on 4>have been dropped for simplicity. 

Premultiplying Equation 3. 2. 4 by $ postmultiplying by 

T- 1 

$ and grouping significant terms produces the following equation 
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(3. 2. 6) 


E 


o 



a 


[h^ E^' 


but 


to = ♦'^t 


(3. 2. 7) 


so that 


b = b"^ $ 
— o — 


(3. 2. 8) 


The first term in brackets is then b and the second is b 

— o — o 

Rewriting Equation 3. 2. 6 making these substitutions yields 


E 


o 




E, 


I 


(3. 2. 9) 


by defining the new variable as 



(3. 2. 10) 


Equation (3. 2. 9) reduces to 




CO b 
— o — o 


T 



(3. 2. 11) 


Equation 3. 2, 11 is of the same form as Equation 3. 2. 1 but the "k” 
subscripts of the latter equation have been replaced by "O" subscripts 
in the epoch form. 
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The"a'in Equation 3. 2. 10 is equivalent to the epoch quantity a^ 
as demonstrated by the derivation. Rewriting the equation for a 

a = Ej^' b + (? (3. 2. 12) 

T 

and substituting for b and b in terms of their epoch quantities, 
and bracketing significant terms yields 

- - (3.2.13) 

t 

The term in brackets is E^ so that Equation 3. 2, 13 reduces to 

a = b ^ E ' b + or^ (3. 2. 14) 

— o o — o 

which is the expression for a^, hence 

a = a (3.2.15) 

o 


3. 3 Explanation of Assumptions Made to Integrate the ESF Equations 

The position and velocity, r (t) and v (t) respectively, in the 
Apollo navigation filter are updated in real time and change continu- 
ously along the path. For the epoch state version of the navigation 
filter the initial position and velocity, r^ and v^, are updated at 
measurement times; however, between measurements these 
vectors remain constant. Similarly, the initial error vector, 
e^, is constant between measurements. Thus between measurements; 
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constant 


(3. 3. 1) 


so the time derivative of the error vector is zero, i. e. 


e 

— o 


= 0 


(3. 3. 2) 


also, since 



(3. 3. 3) 


the time derivate of the epoch covariance matrix is zero 





e e 
— 0—0 


(3. 3. 4) 


Therefore, between measurements 

E = 0 

o 


(3. 3. 5) 


Unlike the current error covariance matrix, E(t), which 
changes continuously along the path, the epoch covariance matrix, 

E^, remains constant between measurements so that it does not 
have to be propagated. In the interval between the (k-1)^^ measure- 
ment and the k^^ measurement, E remains constant or E ' (t, ) = E (t, - 1) 

o o k o k 

as illustrated graphically in Figure 3. 3. 

I 

For ideal two body motion the above results may be applied 
exactly. Also, <^, the transition matrix, may be calculated 
analytically for a conic path. 
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With a disturbing acceleration, a^, present, the concepts 
derived above still hold, i. e. the epoch state and epoch covariance 
matrix may remain constant within measurement intervals. But 
motion under the influence of a disturbing acceleration is not two 
body motion. The equation governing the motion of the spacecraft 
is then 

r + — ^ r = a^(r) (3.3.6) 

r"J 

In order to use the two- body formulation of the navigational problem, 

the disturbing acceleration is considered to perturb r^ and v^ from 

their nominal two body values. Similarly, and thus are 

perturbed from their ideal two body values. However, because these 

change only slowly for non- ideal two body motion, the perturbations 

are ignored and E^ is not propagated between measurements. Proof 

that E^ for the actual path is close to E^ for the conic path and 

varies only slowly is given here, thus validifying the approximation 

that E = 0 for the actual path. Remember that E = 0 for the 
o ^ o 

conic path is exact. 

Let the "c" subscript on the state transition matrix $ and the 
* superscript on the covariance matrix E denote the values of these 
quantities for the conic path. See Figure 3. 4 for an illustration 
of the conic and actual path and their related quantities. 

The differential equation for ^ as given in Section 2. 2 is 


d^ 
d t 


= F$ 
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(3. 3. 7) 




Figure 3. 4 Illustration of Actual Path to Which the Conic Path 
of the Spacecraft is Matched 
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where 


F = 


O I 

G O 


(3. 3, 8) 


For the conic path, the differential equation is 


— = F a. (3.3.9) 

dt ^ 


The current covariance matrix E(t) in terms of the epoch covariance 
matrix for the actual path is given by 


E(t) = $E 

o 


but 


E"(t) = E(t) 


(3. 3. 10) 


(3. 3. 11) 


so 


$ E 
c o 




= $ E $ 
o 


(3. 3. 12) 


Solving for E^' in Equation 3. 3. 12 and grouping significant terms 
yields 


V = (3.3.13) 

Define the new variable as 
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(3. 3. 14) 


$rr-, = ^ $ 

C 


SO that Equation 3. 3. 13 expressed in terms of $,p is 


^o' " ^o 


(3. 3. 15) 


$ expressed in terms of $,j, is 




(3. 3. 16) 


Taking the time derivative of 3, 3. 16 yields 


• • 

Substituting Equation 3, 3. 9 for and 3. 3. 7 and 3. 3, 16 for 
Equation 3, 3, 17 becomes 

Fc $c ^ ^c-^T " ^ (3.3.18) 

Collecting like terms Equation 3. 3. 18 reduces to 


*T 


Solving for <>,j, 


yields 


But 


- F^) «>J 






T 


(3. 3. 19) 


(3. 3. 20) 


(3. 3. 21) 
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(3. 3. 22) 


Ft = % (F- 


Rewritting F - F^ in terms of its constituents 


F - F 


O I 
G O 


O I 


Gc O 


and simplifying replacing (G - G^) by G^^, F - F is given by 


F - F 


o 

O 

^ad 

O 


Defining F - F^ as F^^ so that Equation 3. 3. 22 becomes 


= V ^ad*c 


Equation 3. 3. 20 is then written as 


c ad ^c 


By observing the components of Equation 3. 3. 26 it is seen that 
is small;also since initially <f> (t^, t^) = I, remains close to 1, 
and from Equation 3. 3. 16 it is seen that 




(3. 3. 27) 
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Also, since 



E 

o 


and is close to I 


E s 
o 



(3. 3. 28) 


The assumption that = O between measurements for 
non ideal two- body motion greatly simplifies calculations. The 
epoch state covariance matrix does not have to be propagated 
within measurement intervals. Also, because $ is close to 
the transition matrix may be calculated agebraically, as in the 
case for ideal two-body motion, as opposed to integrating a differ- 
ential equation for 


Analytical calculation of the state transition matrix, $, is 
explained here. Consider the state vector x wher^ 



(3. 3. 29) 


r and v are functions of initial position and velocity 


r = r (r^, v^, t) (3. 3. 30) 

V = J (£q» _Jq. t) (3,3.31) 
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and t is the independent variable upon which and thus r and v 

are dependent. Taking partial derivatives of r and v at time yields 


9 r a r 

6r, = 6r + -^ 6v (3.3.32) 

— k _o ^ — o 

— o — o 



9y: 

^£o 




9 V 
+ _H_ 

9 V 
— o 




(3. 3. 33) 


Expressing Equations 3. 3. 32 and 3. 3. 33 in matrix form results in 



9r 9r 



9 V 9 V 


9 r^ 9 V 
— o — o 


6 r 
— o 

®Io 


(3. 3.34) 


Defining the 6x6 matrix of partial derivatives with respect 
to r and v as #(t, , t ) Equation 3. 3. 34 becomes 

O — O cC O 

6 fk = $(t^, t^) 6x^ (3. 3. 35) 

The partial derivatives of the state transition matrix as given 
in Reference 2 are written here 


9 r 

9 r 
— o 


v(y 


(3. 3. 36) 
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(3. 3. 37) 


9 V j- rp t-p_ -I 

_n. = C(t)V'"(t - R(t)^ 


9 r ^ ,T 

_Z- = R(t) = -R (t 

3v 
— o 


(3. 3. 38) 


where 


9 V 




V(t) 


(3. 3. 39) 


C(t) = 


(3. 3. 40) 


The R and V matrices are given by 


R(t) =— {[U„ (r - r^) + c - U2 (v - v^) r^”^} + G I (3. 3.41) 




V(t) = i- (U„ r r ^ - c r V ^) + — { v - v^) (v - v^)"^ + (3. 3. 42) 


where 


/JIc = 3U5 - X - Ug ^/M(t - t^) 


(3.3.43) 


and X is obtained by solving Kepler's equation. The R and V 
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matrices obtained from Equations 3. 3. 41 and 3. 3. 42 by interchanging 
X and t by -X and -t as well as interchanging r, v and r^, are 

R ' (t) = — {[U 2 (£q - £) - ^£ 0 ^ ” ^2 -o ” - G I (3. 3. 44) 

V*(t)=J— (U„r r'^ + cr^ v'^) + il (v - v) (v - v)”^ + F I (3.3.45) 
^o M 

Because of the assumptions allowing for the analytical 

calculation of computer run time and computer storage space 

are conserved. This approximation is not without some loss of 

accuracy, however. The variation of parameters equations for 

r^ and v^ derived in section 3. 4 are exact but the way they are used 

introduces some error in That is, the position and velocity of 

the spacecraft are matched with the position and velocity of a 

conic path yielding a conic epoch position and velocity that differ 

from the actual r and v . 

— o — o 

3. 4 Derivation of the Variational Equations 


The variational parameters for the epoch state filter are the 
epoch state variables r^ and v^, and the true anomoly difference 0. 
Derivation of the variational equations in terms of these parameters 
is presented here. 

Current position and velocity vectors r and v can be expressed 

in terms of their values r and v at some epoch time t as follows 

— o — o ^ o 
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defining the matrix of scalar quantities F, G, F^, as Y where 


F = 1 - -1 = J_ (r - aU ) 

r r 
o o 


G = r Uj - cr U2 


The U functions used in this derivation are given in terms of 6 as 


ar r 

= 1 2 (1 - cos e ) 

° p 
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r 

U, = sin 6 (1 - cos 6) 

/p P 


r r 

U2 = ? (1 - cos 0) 

P 


r r 

Uo = />/]i (t - t ) — sin 0 

^ ° ./P 

with 


(3. 4. 9) 


(3. 4. 10) 


(3.4. 11) 


r 


1 + - 1) cos 0 - 




sin 0 


(3. 4. 12) 


P 


2 r_ - O' r. 



(3.4. 13) 



(3. 4. 14) 


and 


a 


o 



(3.4. 15) 


Equation 3. 4. 1 may be solved for and by premultiplying both 
sides of the equation by Y ^ thus obtaining 



(3.4. 16) 
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where the determinant of the matrix Y is unity. The perturbation 
derivatives of r^ and v^ may now be calculated from Equation 
3.4. 16 applying the formal rules for variation of the orbital elements ' 
Briefly this means Equation 3. 4. 16 is to be differentiated according 
to the usual rules of differentiation but r is treated as a constant 
and the orbital elements as variables. The term dv/dt is replaced 
by a^ and dx/dt by df/dt where ^ represents the change in x arising 
solely from changes in the orbital elements due to a ^ . 

Formal differentiation of Equation 3. 4. 16 yields 



(3.4. 17) 


(3. 4. 18) 


The upper left-hand element of 3. 4. 18, is obtained by 

d t 

formally differentiating Equation 3. 4.6. Thus, 


dG^ j dU2 (x; a) 

d t r d t 


(3. 4. 19) 
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dU2 (x; O') 

where is derived in Reference 2 using the formal rules 

d t 

for variation of the orbit elements and is written here as 


d t 


U, 


+ (1/2) U, 

d t 


acf 
d t 


-d/2)U, 


2 d a 
‘ dt 


(3. 4. 20) 


The perturbation derivative for a , the reciprocal of the semimajor 
axis, is given by 



(3. 4. 21) 


In order to express Equation 3. 4. 20 in terms of and F^, the 
first term on the right-hand side is multiplied' and divided by 
JJi /r^ and— U 2 Z ’ —d kidded to and subtracted from the second 
term so that 


dG 


t 


d t 


^ (11+ 1 / 2 -^ dg ^ r_ JJi ^1 - 

JJl dt u dt^ror” 


J_ 

nr 



Z -Zd 


+ — U9 V • a . - — U„ V • a . 
fi 2 - -d 2 - -d 


(3. 4. 22) 


The term in brackets is F^, and the third and fourth term simplify to 

~ U„ V • a . G, so Equation 3. 4. 22 becomes 
(,l 4 — —a t 
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d G ^ 

i = J- U„ V • a , G + — (-ii + 1/2 Uo — ) F, (3. 4. 23) 

dt M dt ^ dt ^ 


- — U V 
M 2- 


(3. 4. 23) 


Similarly, differentiating the expressions for G, and F the following 
Equations are obtained 


dG 
d t 


U 


^ V • a^ G - — (11 + 1/2 U3 -1^ ) F 


a/F dt 


d t 


+ — r • a , 

Ki - -d 2 


dF 


^ (A1 + Jl U 3 -1^ ) G - 1- (v - v) a . F 

2 9 ^ ^ 4 . ^ pi O a t 


dt r dt 2 dt 
o 


dF 
d t 


(II + 1/2 Uo — ) G + 1 (v - v).a. F 

r 2 dt dt n 

o 


u 


(v - v) • a , F 
— o — — d 


(3. 4. 24) 


(3. 4. 25) 


(3. 4. 26) 


Details in the derivation of Equations 3. 4. 24 to 3. 4. 26 are given 
in Appendix C. 


Equations 3. 4. 23 to 3. 4. 26 can be expressed in matrix form 
as a matrix multiplied by Y ^ plus another matrix, that is 
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d Y 


-1 


d t 


1 TT T 

3d 




+(l/^Uo-^ 
dt dt 


o 

7^ 


^ +(l/2)Uo 
dt dt 


r , T T. 

- 3 ) 3d 




+ 


1 TT T 

— ^2 3 3d 


O 


i U, a, 

p. 2- -d 


r , T T, 

— (v - V ) a 
^ -o - -d 


( 3 . 


The variational equation for the epoch time, t^, may be 
calculated from Kepler's equation 


./5T (t - t^) = Uj + Ug = rj Uj - cUg + U, 


( 3 . 


by formal differentiation, treating t as a constant and using the 
variational equations for the time derivatives of the U functions 
given in Reference 2. This equation is 


d t 




- r 


d t 


^ +{1/2)V^ — 
d t d t 


2 d t d t 


( 3 . 


where c has been defined by 


= 3 Ug - X U 4 - Ug Vm (t - t^) 


( 3 . 


4 . 29 ) 


4 . 28 ) 


4 . 29 ) 


4 . 30 ) 
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If t varies so that +(l/2)U„-^ = 0, Equation 3.4.29 reduces to 


O' _ 


d t 


d t 


d t 


d t 


U 




2 d 


CT 


d t J]I d t 


(3. 4. 31) 


When— ^ is replaced by Equation 3. 4. 21 and-^^ by the following 
d t d t 

variational equation 


dq _ 1 


£ • 3d 


dt /m 

the variational equation for the epoch time becomes 


(3.4. 32) 


d t 


dt M 


(c V + U 2 £ ) • _a^ 


dY 


-1 


Getting back to • 


,when-^ +(1/2 )Uo-^3. 


(3. 4. 33) 


= 0 is substituted 


d t dt d t 

into Equation 3. 4. 27, this matrix simplifies to 


d y“ ^ 

1 TT T 

-'^21 3d 

0 


^t 

-G 

d t 

0 

r , T Tx 

-<Jo -3 >3d 


-^t 

F 

.i.. 


+ 


- Up v'^a , 
ju 2 _ _d 


O 


1 T 

3 u., r^ a . 

jU 2 


r / T T. . 

-- <Zo -I >^d> 

M’ 


(3. 4. 34) 
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Taking components of Equation 3.4. 17, the perturbation 


derivatives of r and v are expressed as 
— o — o ^ 


d r 
— o 

dt 


^*^t T dG T _ 
r V - Ga . 

dt “ dt “ 


(3. 4. 35) 


d V 
— o 

d t 


dF 


t T , dF T ^ 
r + — V + F a . 

dt ~ dt “ 


(3. 4. 36) 


Multiplying out the matrices of Equation 3. 4. 34 produces 


d Y 


-1 


d t 


1 TT T r- 

■p I G( 


— (v - v^} a . F, 
pi -o _ _d t 


1 T 
- jL u„ V a . G 
u 2 - _d 


^ /tj* ^ 3^ \ XT' 

-p<Io -I 


+ 


1 TT T 


o 


1 TT T 

71^21 


r , T T. 

- — (v - V ) a , 
— o — — d 

pi 


(3. 4. 37) 


d Y ^ 

When the related components of from Equation 3. 4. 37 

d t 

are substituted into Equation 3. 4. 35, the result is 


^-o _^2rT rrT T „T T T 

Iv a , G, r - V a . Gv - v a . r 
— — d t— — — d — — — d — 


d t 




4- T T, „ 

+ r V J - G a 


-d 


(3. 4. 38) 
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T T T 

Collecting like terms and making the substitution r - Gv 


produces the following equation 



d t 




a , r 
— d — o 






G a 


(3. 4. 39) 


Finally, taking the transpose or Equation 3. 4. 39 and collecting 
like terms results in the variational equation for r^ 



(3, 4. 40) 


- 1 


d y 

Similarly, when the related components of from 

d t 

Equation 3. 4. 37 are substituted into Equation 3, 4. 36 the result is 


d V 


— o rr/T T> ■c’Tj./T T, 

= — I - (v - V ) a . F. r + (v - v ) a , F v 

-o — — d t— — o — — d — 


dt 


(Jo ^ ^^d 


(3. 4. 41) 


T T T 

Making the substitution = F^ r + F v in Equation 3. 4. 41 
yields 


d V 

d t fj. 


° = — [(Vq^^ - v"^) ^d Vq'^ - (Vo'^ - v'^) ^d v^] + F^d (3.4,42) 


Transposing Equation 3. 4. 42 results in the following variational 
equation for v . 
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(3.4.43) 


Using the following relation derived in Reference 2 


li +(l/2)U„-i5L = 
dt dt 


a/^^o 

h 


d6 ii/ll I G , ^2 1 fo A aa\ 

- ( — _ h X r + r) • a . (3. 4. 44) 

dt r h"* “ “ a ~ 


where G obtained form Equations 3. 4. 4, 3. 4. 9, and 3. 4. 10 is 

given by 


G = ^ sin 6 

h 


d 5 

and is determined according to 

d t 


d5 ^ _ _h 

dt dt r^ 


(3. 4. 45) 


(3. 4. 46) 


6 being the variation in the true anomoly difference 0. Equation 

3.4.47 is obtained when-^^ +(l/2)Uo-^-^ is set equal to zero and 

dt dt 

substitutions 3. 4. 45 and 3. 4. 46 are made in Equation 3. 4. 44. Thus, 


//jlr 

0 = ^ ( 


U. 


- ^ sin 0 h X r + — - r ) • a , (3. 4. 47) 

h dt r'^ r h-^ 


d 0 

Solving for in the above equation results in 

d t 


d0 _ h 


d t 


+ — [sin0hxr + 


3 

jUrr 


r] • a 


-d 


(3. 4. 48) 
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Recalling Equation 3. 4. 10 for U2 where 


p = — (3.4.49) 

The second term in the brackets is 
3 

hU„ 

r = h (1 - cos 9 ) £ (3. 4. 50) 

/irro 

Therefore, the variational equation for the independent variable 
0 is 

(3. 4. 51) 

Equations 3. 4. 33, 3. 4. 40, 3. 4. 43, and 3. 4. 51 are the variable 
epoch form of the variational equations. This form of the variational 
equations was used in this thesis for simplicity, since the term 

+ — U„ was eliminated by forcing the epoch time, t , 
d t 2 d t , ° 

to vary according to: t = — (c v + U„ r) ' a^. The variational 

O ^ ^ — — Q 

equations for the fixed epoch case are obtained by setting t'^ = 0 in 
Equation 3 4. 29. 

3. 5 Effect of a Measurement on the True Anomoly Difference 0 

When a measurement is taken, current time is essentially 
stopped and the epoch time, t^, remains fixed. Holding the epoch 
time constant during a measurement incorporation in the variable 
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epoch case does not introduce inconsistencies from the time- of- flight 
standpoint, since the state transition matrix relates variations 
in the state at the given to state variations at the current 

time "t". However, since 0 is the angle between the epoch and 
current position vectors and because measurement incorporation 
changes the estimate of these vectors by 6 jq 6^ respectively, 
the true anomoly difference may change by an amount ^0. This 
is illustrated in Figure 3. 5. Tie epoch variation in the true anomoly 
difference is given by 


6 6 , 


ieo' 


(3. 5. 1) 


where_i - is the unit vector in the direction of the epoch change in 0 

and is normal to r^ . Similarly, the current variation in the true 
anomoly difference is given by 


60 = 



r 


(3. 5. 2) 


The total change in 0 is the difference between the current and the 
epoch deviations, that is. 


^6 = 60 - 66^ (3. 5. 3) 

Substituting Equations 3. 5. 1 and 3. 5. 2 into 3. 5. 3 produces 
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(3. 5.4) 


A0 = 



6 r 


r 


A 



Expressing Equation 3. 5.4 in terms of the current and epoch state 
deviation vectors yields 


A0 = 



(3. 5. 5) 


A 

The estimate of the current state deviation vector 6x is determined 
from measurement data by 

= w 6 q (3.5.6) 

Extrapolation of the epoch state deviation vector to the current time 
is accomplished by use of the state transition matrix according to 



6x = <i>{t, t^) 6 Xq 


(3. 5. 7) 


where 6x is determined from measurement data by 
— o •’ 


6x, 




(3. 5. 8) 


Making these substitutions in Equation 3, 5. 5 and factoring out 
6 q results in the following Equation for the total change in the 
true anomoly difference: 
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( 3 , 5 . 9 ) 
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3. 6 Comparison of the Conventional and Epoch Formulations of the Navigational Problem 


CONVENTIONAL FORMULATION 


EPOCH FORMULATION 


STATE VARIABLES AND VARIATIONAL 
EQUATIONS OF MOTION 


r(t), v(t) 

i = If 

V = - [f(q) r + 6] + a^(r) 


osc 


where 


6 • (6 - 2 r) 


q = 


f(q) = 


q (3 + 3q + q ) 


1 + (1 + q) 


3/2 


r (t), 

v^(t) 


— o 

— o 


r = 

- U 

[ (£o ■ 

— L) 

u 


• 

I (v 

- v) (v^ - v) 

— o 

u 

— — o — 




e = \ + -i-[sinehxr+h(l - cos 6) r] • a , 

for which the epoch time, t^, is forced to vary 
according to 


t = — (cv + U„ r) • a 




2 A' -d 


where 


U. 


r r 


(1 - cos 6 ) 




E = 
where 

03 = - 

a = I 


CONVENTIONAL FORMULATION 


EPOCH FORMULATION 


OPTIMUM LINEAR ESTIMATION EQUATIONS 


6x'+o3(6q-b 6x') 


T ^ 

:i - 60 b ) E 


E' b 


E' b + O' ^ 
— n 




(1-03 b E ' 
rro — o o 


E ' b 
o — o 


E 


where 


03 
— o 


T ' 2 

a =* b E b + O' 
o —00—0 n 


b^ E' b + O' ^ 
— — n 


= a 


b = ^ b 
— o — 


< X I 



3. 7 Statistical Equations for Error Using the Epoch Formulation 


of the Navigational Problem 

When the epoch formulation is used rather than the conventional 
formulation of the navigational problem, error is introduced because 
of the assumption that allows for the analytical calculation of the 
state transition matrix as opposed to integrating a differential equation. 
The statistical equations for this error are developed here. Let <i>^ 
denote the conic state transition matrix calculated for the epoch 
state filter. The current estimate of the state deviation vector for 
this filter is given by 


6x = co6q 


(3. 7. 1) 


where 


E' b 

w = ^ (3.7.2) 

a 

and the estimate of the epoch state deviation vector is given by 

6 X = Wq 6 q (3. 7. 3) 


where 





(3. 7. 4) 


The estimate of the current state deviation vector for the epoch state 
filter given in terms of epoch quantities is 


57 





(3. 7. 5) 


6x = 6x^ = 6q 

Let the error vector, incurred using the epoch formulation be the 
difference between the current state deviation vector of the conventional 
filter and that of the epoch state filter. This error is given by 


®d ~ w6q-«i>WQ6q (3.7,6) 

Factoring out 6q and substituting equations 3. 7. 2 and 3. 7. 4 for 
03 and 03 q respectively yields 

Substituting for b^ in terms of b where 

b^ = ^'^b (3.7.8) 

and factoring out b results in 

-Sd " " ^^o' (3.7.9) 


The error covariance matrix, E^, introduced by using the epoch as 
opposed to the conventional formulation of the navigational problem is 



(3. 7. 10) 
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(3. 7. 11) 


T 

[ E' - « E ' b 6q^ [E' - «■ E ' t"'’] 

Ed = 2 2 

The total error vector, e^., the epoch solution of the navigational 

problem is the error inherent in the exact solution of the navigational 

problem, e^^, plus the error incurred using the epoch formulation, 

e .. That is, 

— d 


®t 


e + e . 
— n — d 


(3. 7. 12) 


Solving for the total error covariance noting that^^ and are 
uncorrelated yields 


®t ®t' 


T , 
e e + 
-n— n 


2 


(3. 7. 13) 




(3. 7. 4) 


Since the first term on the right hand side of equation 3. 7. 14 is E^, 
the error covariance matrix inherent in the exact solution of the 
navigational problem, and the second term on the right is just E^, 
then the total error covariance matrix, E^, is given by 


E, = E + E . 
t n d 


(3. 7. 15) 


Thus, the total error in the epoch solution of the navigational problem 
is the sum of the error inherent in the exact solution plus the error 
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introduced by using the epoch as opposed to the conventional or 
exact formulation of the navigational problem. 

The percentage of error introduced in the solution to the navi- 
gational problem by using the epoch state filter rather than the conven- 
tional state filter is obtained by comparing the estimated state vectors 
for both filters with the rms error in the estimate. For a particular 
solution, the percentage of error is given by the following; 


% error 


approximate solution - exact solu tion 
exact solution 


( 3 . 7 . 16 ) 


The magnitude of the error between the position vector determined 
by the epoch state filter and the position vector determined by the 
Apollo navigation filter is given by - £csjrl • A measure 

of the error in the position estimate of the exact or conventional 
solution to the navigational problem is given by the rms estimated 
position error. The percentage of actual error in the position 
estimate introduced by using the epoch rather than the conventional 
state filter is given by 


% actual error in position estimate = ■ TRUE ^ CSF —TRUE ^ 

(rms position error)^^gp 

( 3 . 7 . 17 ) 

Similarly, the percentage of error in the velocity estimate of the 
epoch state filter as compared with the conventional navigation filter 
is determined by 
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% actual error in velocity estimate = 


l-ESF -TRUE^ -CSF -TRUE^^ 
(rms velocity error)^gp 


(3. 7. 18) 

When the percentage errors given by equations 3. 7. 17 and 3. 7. 18 
are small, the ESF may be used in place of the CSF. 
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CHAPTER IV 


COMPUTER SIMULATION RESULTS 


4. 1 Simulation Data 

The purpose for computer simulating the ESF was to determine 
the error in this filtering technique as compared with both the true 
state and the error in the conventional navigation filter. A circular 
earth orbit with a 100 nautical mile altitude and a disturbing accelera- 
tion due to the Jg term of the earth's gravitational potential was 
chosen for study. As a further demonstration of the characteristics 
of this filter, a circular orbit of radius equal to twice the equatorial 
radius of the earth was also studied. Finally, the ESF was simulated 
for the 100 nautical mile orbit with disturbing acceleration due to 
IOJ2 to study the effects of larger disturbing accelerations. 

Measurement data was incorporated at intervals of 10 “around 

the orbit for 80 measurements, thus testing the filter for 800° or 

more than two revolutions. The measurement vector b was alternately 

chosen to be a unit vector in the x, y, and z directions respectively 

for sets of three measurements. This was done so as not to bias the 

problem in any one direction. The error in the measurement, a, was 

2 6 2 

produced by a random number generator with a variance ty of 10 m . 

The initial covariance matrix was chosen to be diagonal with an 

2 

rms position error or 8. 84 x 10 m and an rms velocity error of 8. 65 
m/sec. A diagonal matrix was used so as not to bias the estimation 
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problem, although for this simulation the initial covariance matrix is not 
critical provided, however, that it is large enough. The covariance matrix 
is gradually reduced by measurement incorporation regardless of its 
initial quantity. 

For both the ESF and the CSF, maximum likelihood filtering 

techniques are employed. No approximations are made to extrapolate 

the error covariance matrix for the CSF as were made for the ESF so it 

might seem that the CSF is more accurate. However the CSF has extra- 

T 

polation errors due to the equation E = $E'$ . Since the ESF is an 
approximation to the CSF, a comparison of the performance of these 
filters was made. The results of this study are given for one Monte Carlo 
run rather than the average of many computer runs. 

4. 2 Integration Techniques 

The integration techniques for both filters were compared on the 
IBM 360 model 75 computer. For the integration scheme of the ESF, 
epoch state variables are used to integrate the variational equations 
whereas for the CSF the current state variables are used. Several 
runs of the integration techniques for both filters were made with 
various integration step sizes. The results were compared with an 
exact solution of the equations of motion. This exact solution was 
obtained for a disturbing acceleration of (-|i/20)(r /r ) so that the 

equation of motion reduced to the following two- body equation: 

r 

r +— n -=- = 0 (4. 2. 1) 

20 

Solution to this equation was obtained analytically. 
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Both filters were run for various integration step sizes. The 
error in significant figures between the states for both integration 
schemes and the exact state of the spacecraft was noted after one 
revolution for different integration step sizes. Also noted was the 
computer run time for both filters. Figure 4. 1 illustrates the 
time, in seconds, for the computer runs of the integration techniques 
used in both filters with different integration step sizes for one revo- 
lution. From this figure it is seen that on the IBM 360 computer 
the integration scheme used in the ESF takes longer to run than that 
of the CSF for the same step size. However, reference to Figure 4. 2 
shows that the integration technique used in the ESF is more accurate 
for the same step size than that of the CSF. After one revolution, 
a 5 significant figure error in the estimate of the state as compared 
with the 12 significant figures of the solution was the accuracy chosen 
for the simulation. For this accuracy the errors in the Encke integra- 
tion scheme did not degrade the solution and rectification was not 
required. 

An error of 5 significant figures implied a 1° integration step 
size for the CSF and a 5 step size for the ESF as seen in Figure 4. 2. 
When the step size was eliminated as a parameter from Figures 4. 1 
and 4. 2 a plot of computer run time verses the error in the integration 
technique for both filters was made (Figure 4. 3). The significance 
of this last plot is that when simulated on the IBM 360-75 computer 
the integration techniques of ESF takes longer to run for 
the same integration stap size than its counterpart, however, 
it is more accurate. It was originally thought that 
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the integration scheme for the ESF would be less accurate and more 
time consuming than that of the CSF but this is not the case. In any 
event, both integration schemes could be used for either filter. 


4. 3 Measurement Incorporation for the Simulation 

Actual measurement data, 6q, is incorporated by the filter to 

improve the estimate of the state. This data consists of the true 

T 

measurement value given by b error in the measure- 

ment, a, and is given in terms of its constituents as 

6q = b + Q/ (4.3.1) 


For the CSF, the estimated and true position are determined by Encke's 

A 

Method where? = £oSC -TRUE “ loCS -TRUE’ 

position error vector is given by r - £-pj^-{jE " .^TRUE’ Similarly, 

the true velocity error v - ^-pj^UE by 0 - £-pj^UE ^bat the 

state deviation vector for the CSF is determined from 

6x = 0 ) (6q - b^ 6x' ) (4. 3. 2) 

or 



For the ESF, the epoch state deviation vector is given by 
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(4.3.4) 







) 


or 



4. 4 Measurement Incorporation for Zero Disturbing Acceleration 


In order to demonstrate that the errors in the simulation study 
of the ESF are indeed the errors introduced because of measurement 
incorporation, the ESF was tested for a 100 nautical mile circular 
earth orbit with zero disturbing acceleration. The results of this 
simulation are given in Figures 4. 4 to 4. 11. For this case, the ESF 
and the exact solution to the navigational problem, the CSF, are 
essentially the same. The statement $ = 0 between measurements 
is not an approximation for the ESF with zero disturbing acceleration 
since the motion of the spacecraft in its orbit is two- body. 


In Figure 4. 4 the magnitude of the difference between the esti- 
mated position of the ESF and the CSF, " ^CSF I ’ plotted 

for 80 measurement intervals. This graph shows, that both filters 
are close for the case of zero disturbing acceleration except for a 
slight random error which grows with time. At the eightieth mea- 
surement this error is only 1. 46 meters. Figure 4. 5 illustrates the 
magnitude of the difference between the estimated velocity of the 

ESF and the CSF, |v™_, - I for 80 measurement intervals. 

’ '-ESF -CSF' 
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Figure 4. 6 illustrates the fact that for zero disturbing acceleration 

both filters are comparable except for a small error. This random 

error has a maximum value of 1. 44 meters for the case of Figure 4. 6 

and value of .55 meters after the eightieth measurement. Similarly, 

reference to Figure 4. 7, a plot of the difference between the actual 

velocity error of both the ESF and the CSF ” 

Iv^^.,-, - for 80 measurement intervals, shows that both filters 

'—CSF —TRUE' 

are comparable except for a slight random error which is . 00182 
m/sec at maximum. This is also the value of error after the eightieth 
measurement. To find the error in the ESF that is in excess of the 
error in the CSF at the eightieth measurement, the difference between 
the actual position vectors of the ESF and the CSF, . 55 meters, 
is compared with the expected rms position error of the exact or CSF 
solution, 523 meters as seen in Figure 4. 8, The percentage of error 
in the position estimate using the epoch formulation of the navigational 
problem instead of the conventional formulation is approximately 
. 11% for the case of zero disturbing acceleration. Likewise, the 
error between the actual velocity errors of the ESF and the CSF, 

. 00182 m/sec, is compared with the expected velocity error of the 
CSF solution, . 55 m/sec as seen in Figure 4. 9. The result is that 
a . 33% error exists in the velocity estimate of the epoch state filter 
in excess of the error in the velocity estimate of the CSF. 

In Figures 4. 10 and 4. 11 the expected rms position error and the 

rms velocity error at the epoch are given for the case of a 100 nautical 

mile circular earth orbit with zero disturbing acceleration. The rms 

2 

position error starts out at its initial value, 8. 84 x 10 meters, and 
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then decreases steadily to a value of 445 meters after 80 measurement 
intervals. Similarly, the rms velocity error decreases from its 
initial value of 8. 65 m/sec to . 52 m sec with measurement incorpora- 
tion. Reference to Figures 4. 10 and 4. 11 shows that the rms velocity 
error at the epoch decreases faster with measurement incroporation 
than the rms position error of the epoch. The fact that these last 
two graphs are for the epoch covariance matrix is made clear by 
noting that between measurement intervals, the rms position errors 
are constant. This is in agreement with the fact that = 0 between 
measurements . 

Between measurements, the errors in the current estimate of 
the state of the spacecraft grow with time. These errors are then 
reduced with the incorporation of new data as is seen in Figures 4. 8 
and 4, 9. However, the errors in the estimate of the epoch state 
can only be reduced. They do not grow with time but are constant 
within measurement intervals, and the incorporation of new information 
acts to only decrease the error in the estimate of the epoch state 
(Figures 4. 10 and 4. 11) 

4. 5 Disturbing Acceleration Due to Jg Term 

The ESF was tested for a 100 nautical mile circular earth orbit 
with a disturbing acceleration due to the J 2 term. The results of 
this test were compared with similar results for the CSF. For this 
case, the magnitude of the estimated position deviation vector, 

I 6 r I , is given for both filters after 80 measurement intervals in 
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Figure 4. 12. Reference to this figure shows that the magnitude of 
the position deviation is approximately the same for the ESF and 
the CSF until the thirty- seventh measurement and after that differs 
only slightly. The same is true of the velocity deviation (6v( 

(Figure 4. 13). The differences exhibited in these figures may be 
attributed to the random type of measurement data for the simulation. 
When the filter has been operating for a while and reducing the error 
in the estimate with measurement incorporation, the estimate of the 
state is more accurate. Since the estimate of the state is more 

I 

reliable, the required state deviation is lessened. This is seen 
in Figures 4. 12 and 4. 13. 

As time increases, the difference between the position deviation 
vectors for the ESF and the CSF grows (Figure 4. 14) and is due to 
the difference between the filter gains. However, this difference is 
small, having a maximum value of about 48 meters and a value of 
about 2 meters at the eightieth measurement interval. Similarly, 
the difference between the velocity deviation vectors of the ESF and 
the CSF (Figure 4. 15) has a maximum value of about . 06 m/sec 
and a value of about . 003 m/sec after the eightieth measurement 
interval. 

Figures 4. 12 to 4. 15 are plots of the state deviation vector and 
illustrate the effect of each measurement incorporation on the filter. 
These figures show what is to be added to the estimate of the state 
because of the incorporation of new data. The accumulated effect 
of measurement incorporation on the ESF as compared with the 
CSF is given in Figures 4. 16 and 4. 17 for the estimated position 
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and velocity respectively. Again, the irregularity in these plots 
is due to the random type of measurement data. The significance 
of Figure 4. 16 is that the estimated position difference between the 
ESF and the CSF remains nearly zero for ten measurements and then 
grows to about 40 meters after the eightieth measurement. Similarly, 
the estimated velocity difference (Figure 4. 17) remains nearly zero 
for ten measurements and grows to .065 m/sec after the eightieth 
measurement. 

The errors given in Figures 4. 16 and 4. 17 are for the difference 

between the estimates of the ESF and the CSF with no indication of 

the true state. At times, the estimate of the state as given by the 

epoch formulation of the navigational problem may be more correctly 

aligned with the true state than that of the CSF. This is because the 

ESF does not have the extrapolation errors of the CSF due to the 

T 

equation E = <I> E^ ^ . Also, the approximation that E^ = $ = 0 
between measurements for the ESF is nearly exact. 


A more significant test of the epoch formulation of the navigational 
problem is obtained by comparing the actual errors for the ESF and 
the CSF. The actual error in the ESF is the difference between the 
estimated state for this filter and the true state, |£pigjp " £tRUE^ ' 
and is given in Figure 4. 18. This is easily obtained since for the 
simulation the true orbit of the spacecraft is known. The actual 
error in the position estimate of the ESF has a maximum value of 
about 1800 meters and reduces to 523 meters after eightly measure- 
ment incorporations. The same is true for the actual position error 
of the CSF, seen in Figure 4. 19. When the 
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results of Figure 4. 19 are subtracted from those of Figure 4. 18, the 
difference between the actual position error of the ESF and the CSF 
results (Figure 4. 20). This last graph is probably the most signifi- 
cant of the simulation. What it implies is that both filters are 
comparable. When the actual error difference is positive, the ESF 
has more error in the estimate of the position vector and the CSF 
is the more correct. Conversely, when the actual error difference 
is negative the ESF is more correct. After the eightieth measurement, 
the actual error difference has a magnitude of 16 meters. As compared 
with the corresponding rms estimated position error which has a 
value of 523 meters as seen in Figure 4. 2l, the percentage of error 
introduced by using the epoch formulation of the navigational problem 
is approximately 3, 1%, This is the percentage of error due to the 
ESF in excess of the error inherent in the conventional solution to 
the navigational problem and is small for all practical purposes. 

Figures 4, 22 and 4. 23 are graphs of the actual velocity error for the 
ESF and the CSF respectively. As seen in these figures the actual 
error in the estimation of the velocity has a maximum value of 
about 3 m/sec and is reduced to about . 2 m/sec after the eightieth 
measurement. More significantly the difference between Figures 
4. 22 and 4. 23 as given by Figure 4. 24 shows again that both filters 
are comparable. The ESF is more accurate than the CSF and vice 
versa. Initially, the actual error difference in the velocity estimates 
for both filters has a value of zero. After the eightieth measurement 
this error increases somewhat randomly to a magnitude of . 048 m/sec. 
When this last value is compared with the corresponding rms velocity 
error which the percentage of error using the ESF to estimate the 
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velocity vector is approximately 8, 7%. Again for all practical purposes 
this error is small. 

One of the reasons that the ESF is an accurate estimator of the 
current state is evidenced by the change in the epoch state from its 
true initial value for eighty measurement intervals. In Figure 4. 26 
the epoch position change, ! (t) - seen to be 390 meters 

at maximum and 170 meters after the eightieth measurement. The 
maximum value, | v^(t) - ^ change after 

the eightieth measurement is . 55 m/sec as seen in Figure 4. 27. 

Both the true epoch position and velocity change are small enough 
to insure the accuracy of this formulation of the navigational problem. 

To insure that the errors in the simulation were not due to the 
integration technique, the osculating orbit was rectified every 180° 
for the 100 nautical mile circular earth orbit with a disturbing ac- 
celeration due to the J 2 term. As illustrated in Figures 4. 28 to 4. 31, 
rectification produced no observable change in the estimate of the 
state since the errors in the integration technique are small enough 
so as not to degrade the solution. 

Results of another Monte Carlo run for the 100 nautical mile 
circular earth orbit with a disturbing acceleration due to the J 2 
term are given in Figures 4. 32 to 4. 47. These graphs especially 
Figures 4. 44 and 4. 47 confirm the result of the previous Monte 
Carlo run, in particular that the ESF and the CSF are comparable 
for this orbit and disturbing acceleration. 

The error in the approximation that ^ = 0 between measurements 
for the ESF decreases with the disturbing acceleration which is a 
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function of l/r^_ Because of this, the accuracy of the filter increases 
as the spacecraft travels farther into space away from the disturbing 
influence of the earth. To illustrate this characteristic of the filter, 
the ESF was simulated for a circular earth orbit with a radius equal 
to twice the equatorial radius of the earth. The results of this 
simulation are given in Figures 4. 48 to 4. 53. 

As seen in Figure 4. 48, the difference between the position esti- 
mates of the ESF and the CSF is considerably less for this orbit 
than for the 100 nautical mile orbit. After the eightieth measurement, 
the position difference for this larger orbit has a value of 11 meters 
as compared with 40 meters for the 100 nautical mile orbit. The 
difference between the actual position error of both filters (Figure 
4. 49) varies randomly having a magnitude of only 4. 2 meters after 
the eightieth measurement as compared with 16 meters for the 100 
nautical mile orbit. The rms position error for the circular orbit 
of radius r = 2r„ (Figure 4. 50) has a value of 500 meters after the 
eightieth measurement. Comparing the difference between the actual 
position errors with this value results in a . 84% error. This is the 
extra percentage of error introduced by using the ESF to extimate 
the state of the spacecraft. A . 84% error is a considerable reduction 
when compared with the 3. 4% error for the 100 nautical mile orbit. 
After the eightieth measurement, the difference between the velocity 
estimates of the ESF and CSF (Figure 4. 51) has a value of . 0048 
m/sec for the orbit of radius r = 2 r^ as compared with . 048 m/sec 
for the 100 nautical mile orbit. The difference between the actual 
velocity errors (Figure 4. 52) varies randomly having a magnitude of 
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. 0048 m/sec for the 100 mile orbit. When the difference between 
the actual velocity errors after the eightieth measurement is compared 
with the corresponding rms velocity error which is . 2 ml sec as seen 
in Figure 4. 53, the percentage of error introduced by using the ESF 
for the circular earth orbit of radius r = 2 r^ is 2. 4%. This percentage 
of error is considerably less than the 8% error in the velocity estimate 
for the 100 nautical mile orbit. These results demonstrate that the 
accuracy of the ESF increases for higher orbits becuase the relative 
accuracy of the approximations increase. 

4. 6 Disturbing Acceleration Due to 10 J 2 

It would be interesting to apply the ESF to the re-entry navigational 
problem. However, for this problem the spacecraft is subject to large 
values of disturbing acceleration. To see if the epoch formulation 
of the navigation problem works properly for disturbing accelerations 
due to terms larger than J 2 , the ESF was simulated for a 100 nautical 
mile circular earth orbit with a disturbing acceleration due to 10 Jg. 
Results of this simulation are presented in Figures 4. 54 to 4. 69. 

In the first of these graphs. Figure 4. 54, the magnitude of the 
position deviation for the ESF and CSF is seen to be the same for 
both filters until the eighth measurement. Although the difference 
there if only slight, it becomes more pronounced as time goes on 
and the error in the approximation for the ESF increases. This same 
result is seen for the magnitude of the velocity deviation (Figure 
4. 55). The difference between the position deviation for the ESF and 
CSF is given in Figure 4. 56. In this case, ! S - 6 ?qsF ' ^ 



maximum value of 490 meters as compared with 48 meters for the 
case of a disturbing acceleration due to and a value of 30 meters 
after the eightieth measurement as compared with 2 meters for the 
Jg case. The velocity deviation (Figure 4. 5 7) has a maximum value 
of . 58 m/sec for the 10 J2 case as compared with . 06 m^sec for the 
J2 case and a value of . 50 m/sec after the eightieth measurement 
for the 10 J2 case as compared with , 003 m 'sec for the J2 case. 

This considerable difference is due to the error in the simplifications 
of the ESF which is large for the case of a greater disturbing ac- 

f 

celeration. The difference becomes more pronounced as time goes 
on and the errors of the filter diverge. 


Similarly, the difference between the position estimates for 
the ESF and the CSF (Figure 4. 58) has a maximum value of 400 
meters for the 10 Jg case as compared with 43 meters for the J2 
case and a value after the eightieth measurement of 275 meters 
for the 10 J2 case as compared with 40 meters for the J2 case. 

The difference between the velocity estimates, . 

(Figure 4. 59) has a maximum value of . 51 m^sec for the 10 J2 
case as compared with . 065 m/sec for J2 case and a value of . 35 
m/sec after the eightieth measurement as compared with . 050 m/sec 
for the J2 case. 


The actual error for the ESF with a disturbing acceleration due 
to 10 J2 as seen in Figure 4. 60 has a maximum value of about 1770 
meters and a value of 470 meters after the eightieth measurement. 
The actual error for the CSF (Figure 4. 61) has a maximum value of 
about 1770 meters and a value of 650 meters after the eightieth 
measurement. 
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Again, the difference between the actual errors for both the 
ESF and the CSF varies randomly (Figure 4. 62) implying that the 
ESF and the CSF are comparable. That is, the ESF is more accurate 
in estimating the spacecraft's position vector just as often as the CSF. 
The difference here is that after the eightieth measurement, the 
difference between the actual position errors has a value of about 180 
meters. When this is compared with 500 meters for the rms position 
error (Figure 4. 63), a 36% error is determined. This is considerably 
larger than the 3. 2% error for the case of a disturbing acceleration 
due to the simple J2 term. However, not only is the ESF working 
for the case of a larger disturbing acceleration but it is comparable 
to the CSF in estimating the position of the spacecraft. 

The same is true for the percentage of error in the velocity 
estimate for the ESF. Although there is a 60% error in the velocity 
estimate of the ESF, the ESF and the CSF are comparable as seen 
in Figure 4. 66. 

One of the reasons that the ESF approximately works for the 
case of a disturbing acceleration due to 10 Jg is evidenced by refer- 
ence to Figures 4. 68 and 4. 69. The magnitude of the difference 
between the estimated and true epoch position vectors |rQ(t) - 
has a maximum value of 400 meters and a value of 250 meters 
after the eightieth measurement. Reference to Figure 4. 69 shows 
that the magnitude of the difference between the estimated and true 
epoch velocity vectors has a maximum value of 2. 5 m 'sec and a 
value of . 55 m/sec after the eightieth measurement. 
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The ESF was simulated for the case of a 100 nautical mile orbit 
with a disturbing acceleration due to 100 J2. However, for this case, 
the ESF worked properly for a little over a fourth of a revolution. 

After that the errors in the approximation that $ = 0 for the ESF 
significantly degraded the solution and the ESF did not work. 

Included in this thesis are the computer programs written in the 
MAC language for the simulation of the ESF and the CSF. 

4. 7 Computation Time on AGC 

The relative execution time for the various operations differs from 
computer to computer. The run time solutions given previously 
for the integration techniques of both filters on the IBM 360 computer 
are not necessarily the same for the Apollo computer. To compare 
the computation time for the epoch and the conventional state filters 
on a spacecraft computer such as the Apollo Guidance Computer, 

(AGC) the explicit computational algorithms for both solutions to the 
navigational problem are given here. However, only an approximation 
to the actual computation time for each of these solutions is deter- 
mined according to the total number of various arithmetic and branch- 
ing operations. Using the information presented here, a comparison 
of both filters for computers other than the AGC may be easily made. 

The equations used in the computer subroutines for the ESF and 
the CSF were described in previous sections of this thesis. In the 
following paragraphs, the sequence of computations for both solutions 
of one navigation cycle are given precisely: 
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EPOCH STATE FILTER 


Input 


b : 
6q: 




~~2 
a : 



measurement vector 
measurement data 
initial epoch time 
initial epoch position 
initial epoch velocity 
initial true anomaly difference 
initial epoch covariance matrix 
earth gravitational constant 
measurement variance 
disturbing acceleration term 


Initialization of Loop 


r 

— o 



^/M = SQRT (fi) 
i = 1 


Iterative Loop 


1 ) 

2 ) 


ie 



i 

— z 


X I 


X 1 
— r 
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3) A6 = — (i-- 6r)--l (i • 6r ) 

r -0 - -0Q -o 

4) 0 = e + A0 

5) r = r + 6 r 

6) Set = 1 

7) Set k = 0 

8) Call DIFEQ to integrate 0, r^, with the following iterative loop 
8a) Cos 0 = cos (0) 

8b) Sin 0 = sin (0) 


8c) r^ 
8d) 

8e) Oq 


8f) a = 

8g) p = 
8h) h = 
8i) r = 


8j) 

8k) Ug 
81) F 


llo' 




r • V 
— o — o 


^o ^ 

2r - a r ^ - a ^ 
o o 


^Aup 


n ^0" 

1 + (i- — 1) COS 9 - — ~ sin 9 


7^1 




^ sin 9 - — - (1 - cos 0) 


^/p 


r r 


— (1 - COS 6) 


U. 


1 - 
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8m) = 

- u, 

r r 

o 

8n) 

G = 

^2 

8o) 

^t = 

1 

r 

8p) 

r = 

F r + G V 
— o — o 

8q) 

r = 

l-l 

r 

8r) 

i = 
— r 

r 

8s) 

cos 0 

= i • i 
— r — z 

8t) 

V = 

F, r + G, V 
t— o t — o 

8u) 

h = 

r X V 


r_ 2 r, 

8 v) = -1. 5-i|- Jg (— ) [ {1 - 5(cos 0 ) } + 2 cos 0 

r 

8w) = -^ + — [ sin 6hxr + h(l - cos 6) r] ' a . 

d 1C 

8x) = J_ U„ [ (r - r) + V r"^] a , - G a . 

dt ^ ^ -o 


8y> — ° = [— <io ■ '3o ■ »d ^3d 

cit jU 


8z) Set k = k + 1 


8Q) If k < 4, Go to 8a 


9) Set ^ i + 1 


10) If < 3, Go to 8a 
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11) 

1 /a 

1 



O' 

12) 

r = 
o 

\r 1 

o ' 

13 ) 

Solve Kepler's equation for x 

14 ) 

^3 = 

^ (1/0?) [x - uj 

15 ) 

^/^■(t 

- - ^0 ^1 + ^^0^2 + U3 



2 

16 ) 

^4 = 

(i/c,)[JL_ - Uj] 



3 

17 ) 

^5 - 

' (1/0) [ ^ - U„] 

6 

18 ) 

c = 

-^ [ 3 Ug - X - U2 (t-T)] 

19 ) 

R(t) 

= _i {[ U2 (r - r^) + c v] - [ Ug (v - 

M 

20) 

V(t) 

= ^ 3 [^^ 2 £)i:o - <^£>£0 ^ ^ ° <£-£0 

21) 


) = — „ [ (U2 r ) r^ + (c r ) v^] + [-£ (v 
To'" M 

22) 

r'^-i 

= (r"^)"^ 

23 ) 

C(t) 

T- 1 T 
= R W 

24 ) 

9 r 

9 r 
— o 

= V (t^)T 

25 ) 

9 V 
S£o 

>i< rp rp 4 

= C (t) V''(t^)^ - R(t)^"^ 




26 ) 


‘o> = ( 



\ C(t) V*(y'^ - R(t)'^"‘ 


R(t) 


V(t) 
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27) 

1 O' 

0 

II 

b 

28) 

a = 
0 

£0 • 

29) 

II 

0 

31 

— E ' b 
a 0-0 

0 

30) 

E = 
0 

E ' - 0 )^ b^'^ E^' 
0 -0 —0 0 

31) 

o> 

1 ><!> 
0 

II 

“0 

32) 

Divide 

6 X into 6 r and 6 v components 
0 —0 —0 ^ 


ie, 6r^ = 6 x^ .... 6 v^ = 6 X 3 ... 

33) 6 x = <E> 6 x^ 

34) Divide 6 x into />r and 6 v components 

35) Set i = i + 1 

36) Go to 1 

CURRENT STATE FILTER 

measurement vector 
measurement data 
initial time 
initial position 
initial velocity 
initial covariance matrix 


Input 
b : 
6q: 

t: 

r ; 

— o 

V : 

— o 
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M: earth gravitational constant 


2 : measurement variance 

a 

J 2 ' disturbing acceleration term 


Initialization of Loop 


r = 
o 



V 

o 



a = 




I fa - — 
a 


a = — (r • V ) 
o /— — o — o 


JH = SQRT (m) 


Set upper left-hand 3x3 corner and lower right-hand 3 
elements of F equal to zero 

Set i = 1 


Iterative Loop : 

1) Set off diagonal terms of $ equal to zero 

2) Set diagonal terms of <i> equal to one 

3) ^ ^ 

4) p = V + 6v 

5) set 1=1 


3 corner 
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6) 

7) 

7a) 

7b) 

7c) 

7d) 

7e) 

7f) 


Set k = 0 

Call DIFEQ to integrate 6 , u with the following iterative loop 
Determine x by solving Kepler's equation 
U 2 (x; O') = (1/a) [ 1 - Uq(x; a)] 

^ = ^o ^o ^ ^o ^1 ^2 


U. 


F = 1 - 


F = — ^ U 
^t 1 

r r_ 


7g) 

Gt = 

1 

2 



r 

7h) 

— osc 

= 

F r + 
—0 

7i) 

V ^ 

—osc 


F, r 
t —0 

7j) 

r 

osc 

= 

osc ^ 

7k) 

r = 

r „ + 6 
-osc — 

71) 

r = 

|r 

1 



r 


7m) 

i = 

K» 




JL 

r 


7n) 

COS 0 

= 

i • i 
— r — 



(6 

- 2r) ' 

7o) 

n = 



4 


v,2 


84 



7p) 

7q) 

7r) 

7s) 

7t) 

7u) 

7v) 

7w) 

7x) 

8 ) 

9) 

10) 

11 ) 

12) 

13) 

14) 

15) 


f(q) = q- 


3 + 3Q+ q 
1 + (1 + q)^'^ 




1 5 M j 

l.t3-2 Jg 
r 




1-5 (cos 0 ) } + 2 cos 


G(t) = JL [(3 r)r'^ - l] 

Set lower left-hand 3x3 corner elements of F equal to their 
respective elements of G(t) i, e. Fj^g = G^, F^^g = G^ . . . 

d6 


d t 
dv 
d t 


= V 


osc 


(f (q) r + 6 ) + a 


-d 


At. = F $ 


dt 

Set k = k + 1 


If k < 4, Go to 7a 


Set t = I + 1 


If -t < 11, Go to 7a 
E = «>E' 


a = b • (E' b) + ex' 


w = - E' b 
~ a ~ 

E = E' - ct' b'^ E' 


6 X = ^ 6q 

Divide 6x into 6r and 6 v components i. e. , 
6r = 6^ 6v^ = 6x„... 
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16) Set i = i + 1 


17) Go to 1 

Ite rative Solution of Kepler's Equation for x 
KEl Set j = 0 

KE2 k jl- ° [3„ 2_ r 

r 2 r br 

o o o 

[>/^T(t - T )] ^ + ... I 

2 , 2>2 

CyX (CVX ) 

KE3 U (x ; o) = [ 1 - — + — - . . . ] 

° ^ 2! 4! 

2 . 2,2 

cyx (o-x ) 

KE4 U, (x^; cy) = [ 1 2- + ...] 

" ”31 5! 


KE 5 

U3 (Xn; O') = (l/o')[x^- “)] 


KE6 

./^(tn-T) = U^(\: a) * cr^ 


KE 7 

■■n “ ^0 Uo<V "0 + 

U2 

KE8 

n/M t - t 

X , . = X - 

n+ 1 n 



X 

n 


KE 9 

X = X . . 

n n+ 1 


KEIO 

j = j + 1 



KEll If 3 < 4 , Go to KE 3 

These equations are given in Reference 3 
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DIFEQ (Common to both techniques) 

Given the differential equation; 
dy/dt = f(t, y) 

where y is the dependent variable, t is the independent variable and 
At is the increment, the value ofyatt = t + At can be obtained from 
the following process: 

DQl Set y = Yq and t = t^ (i. e. , their initial values) 

DQ2 Set At = h 

DQ3 Evaluate dy/dt = f (t^, y^) 

DQ4 Evaluate k, = hf(t , y ) 

1 o’ ‘’o 

DQ5 Set t = t^ + (h/2) and y = y^ + (kj/2) 

DQ6 Evaluate dy/dt = ^ (t^ (t»/2), y^ + (k^/2)) 

DQ7 Evaluate k„ = hf(t + (h/2), y + (k, /2)) 

DQ8 Set y = y^ + (k2/2) 

DQ9 Evaluate dy/dt = f (x + (h/2), t + (k„/2)) 

O O u 

DQIO Evaluate k„ = h f (t + (h/2), y + (k.,/2)) 

DQll Set t = t + h and y = y + k„ 

DQl 2 Evaluate dy/dt = f (t^ + h, y^ + k^) 

DQl3 Evaluate k. = hf(t + h, y + k„) 

DQ14 Evaluate k = (k^+2k2+2k3 + k^)/6 
DQl 5 Set y = y^ + k 
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y at this point contains the desired result 
These equations are given in Reference 9. 

Since the time consumed for the input and initialization of 
the loop is a small part of the total computation time for both algorithms, 
these parts are ignored in the following calculations. Table I contains 
a list of some operations common to both filter algorithms and the 
number of functions involved in these operations. The number of 
various arithmetic and branching operations required by each line 
of the iteration loops for the ESF and the CSF algorithms are given 
in Tables II and III. The parameter "i" represents the number of 
iterations performed over the whole iterative loop. Parameters 
"-t" and "k” represent iterations performed over the extrapolation 
of the state. The maximum value of "l" is 2 for the ESF and 10 for 
the CSF and the maximum "value of "k" is 4 for both filters, "j" 
represents iterations performed over the solution to Kepler's 
equation and has a maximum value of 4 for both filters. 

It is pointed out that in formulating the algorithms, little 
attempt has been made to organize the computation so as to minimize 
the overall execution time on the AGC. Tables II and III and sub- 
sequent tables derived from these, represent a reasonable count of 
the number of various operations required to execute the algorithms 
for both filters. In Tables II and III computation has been divided 
into three sections: (A) extrapolation of the state, (B) extrapolation 
of the covariance matrix, and (C) update of the estimates. The 
computation time for the first and third of these sections should be 
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approximately the same for both filters. However, if the method of 
extrapolating the state for the ESF is too time consuming, then 
this method can be replaced by that used in the CSF with a slight 
modification. The big savings for the ESF comes in the extrapolation 
of the covariance matrix. Isolation of the computation time for this 
section emphasizes this savings. 

The total number of various arithmetic and branching operations 
required for each of the algorithmic divisions described previously 
i. e. , (A) extrapolation of the state, (B) extrapolation of the covariance 
matrix, and (C) update of the estimate, are presented in Tables IV 
and V for the ESF and the CSF respectively, after one complete 
execution. 

The relative time required for a single execution of each of 
these operations on the AGKD execution of each of these operations on 
the AGC is summarixed in Table VI. The information is adopted 
from the Users Guide to the Block II AGC/LGC Interpreter, 

Reference 11. 

The total number of "add times" required by one complete 
execution of the ESF and CSF algorithms is given in Table VII in 
"i" iterations. Finally a rough estimate of the overall computation 
time for the AGC required by the ESF and CSF algorithms is given 
in Table VIII for the algorithmic divisions previously described 
and their total. These tables give the results of converting the infor- 
mation for "our case" in Tables IV and V into actual computation 
time in seconds. The results given in Table VIII were expected. 
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and in particular those of Section B. Even if the integration scheme 
for the ESF, Section A, were more time consuming than that of the 
CSF, the latter method could be substituted in its place. However, 
there is a considerable savings in the execution time for the ESF 
because of the assumption that $ = 0 between measurements. This 
savings is evidenced by Section B of Table VIII. The assumption 
allows for $ to be calculated analytically rather than by integrating 
a differential equation. When $ is calculated by integrating a differen- 

9 

tial equation using the MAC subroutine DIFEQ , the integration is 
carried out in four steps to yield <|;for 1° increments. All of the 
equations used in determining the differential equation for 4 are 
sequenced on four times to determine $for the 1° increment. To 
calculate $ at 10° measurement intervals, the equations for 4 are 
cycled forty times. Using the assumption E = 0 between measurements 
for the ESF, 4 > is calculated only once for the 10° measurement 
intervals and need not be computed along the path. 
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Operation Add Multiply Divide SQRT Initialize 


(3x3) Matrix Transpose 





6 

(3x3) Matrix Inverse 

14 

27 

9 


6 

(3x3) (3x3) Matrix Matrix Multiplication 

18 

27 




(3x3) (3x1) Matrix Vector Multiplication 

6 

9 




(3xD (1x3) Vector Vector Multiplication 


9 




(6x6) (6x6) Matrix Matrix Multiplication 

180 

216 




(6x6) (6x1) Matrix Vector Multiplication 

30 

36 




(6x1) (1x6) Vector Matrix Multiplication 


36 




(3x1) Vector Cross Product 

3 

6 




(3x1) Vector Dot Product 

2 

3 




(6xD Vector Dot Product 

5 

6 




Magnitude (3x3) Vector 

2 

3 


1 



The transpose of a (3x3) matrix is considered as 6 initializations since only the off-diagonal 
elements change locations 

TABLE I 

Functions Used in Some Operations Common to Both Filter Algorithms 












Line 


Section 


Add 

Subtract 


1 

c 

2 

c 

3 

c 

4 

c 

5 

c 

6 

A 

mm 

A 

8 

A 

9 

A 

10 

A 

Di 

B 

12 

B 


B 1 

1 


14 


B 


li 



9i 


4 


5i 




























































































































CD 



(Contd. ) 



Trans. 

If 

? 

SQRT 

Function 

(Branch) 

Initialize | 




























































TABLE II (Contd. ) 


■ 

Section 

Add 

Subtract 

Multiply 

Divide 


DQ7 

A 


7/{,i 



DQ8 

A 

7^i 


Ui 


DQ9 

A 

77^i 

149^i 

22li 


DQIO 

A 


Ui 



DQll 

A 

8^i 






77-li 

149-ti 


1 

DQ13 

A 


7-ti 

1 ' '■ ■■ ' < 

1 

DQ14 

A 

2U 

I 7^i 

1 


DQ15 

A 

7-ti 

1 


DIFEQ loop used to evaluate DQ3, 

8a 

A 





b 

A 





c 

A 

2kli 

3k-ti 



d 

A 

2kli 

3kli 



e 

— 1 

A 

2k<ti 

3k-ti 

Ikti 

— 



Ik-ti 


Ik-ti 


2kli 




Trans. 

If 

; 

SQRT 

Function 

(Branch) 

Initialize j 

















































































































































TABLE II completed 


Line 

Section 

Add 

Subtract 

Multiply 

Divide 

SQRT 

Trans. 

Function 

If 

(Branch) 

Initialize 

w 

A 

lOkti 

18k^i 

2k^i 





X 

A 

I6k-ti 

36k^i 

Ik-ti 





y 

A 

15k-ti 

24kti 

Ik-ti 





z 

A 

lk-t,i 







Cl 

A 






Ik-ti 


Solution to Kepler’s equation used to determine 13 

KEl 

B 







li 

KE2 

B 

li 

2i 

li 


li 



KE3 

B 





Iji 



KE4 

B 





IJi 



KE5 

B 


Iji 






KE6 

B 

2ji 

2ji 






KE7 

B 

2ji 

2ji 






KE8 

B 

2ji 

2ji 

Iji 





KE9 

B 







Iji 

KEIO 

B 

IJi 







KEll 

B 






Iji 



































Line 


Section 


Add 

Subtract 


Multiply 




CD 

00 






1 

C 



2 

C 



3 

C 



4 

C 

3i 


5 

A 



6 

A 



m 

A 

3,400i 

4, I60i 

B 

B 

10, 080i 

11, 480i 

8 

B 

l^i 


9 

B 



10 

B 

360i 

432i 

n 

B 

36i 

42i 

12 

B 

30i 

36i 

13 

B 

66i 

72i 

14 

C 


6i 

15 

C 




TABLE III (CSF) 




































































TABLE III (CSF) (Contd. ) 


Line 

Section 

Add 

Subtract 

Multiply 

Divide 

— 

SQRT 

Trans . 
Function 

If 

(Branch) 

Initialize 


16 

C 

li 







17 

C 






li 



"diFEQ for 6, V and 6 

DQl 

A 







7-ti 


B 







36-ti 

DQ2 

A 







Hi 

DQ3 

A 

74^i 

98<ti 

17-ti 

3li 


5-ti 

8U 


B 

189t.i 

251^i 

l^i 




9li 

DQ4 

A 


6-ti 








B 


36-ti 







DQ5 

A 

7^i 


6t,i 







B 

36-ti 


36/1 






DQ6 

A 

74-ti 

98-M 

17-ti 

3-ti 


5-ti 

8-a 



B 

189fi 

251^i 

1-ti 




9-ti 

. 
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TABLE III (Contd. ) 


Line 

Section 

Add 

Subtract 

Multiply 

Divide 

SQRT 

Trans. 

Function 

If 

(Branch) 

Initialize 

DQ7 

A 


6-ti 







B 


36-ti 






DQ8 

A 

6-ti 


6-ti 






B 

36^i 


36-ti 





DQ9 

A 

74-ti 

98^i 

17^i 

3^i 


5-ti 

8-ti 


B 

189ti 

25l>ti 

Hi 




9-ti 

DQIO 

A 


6-ti 







B 


36-ti 






DQll 

A 

7U 








B 

sen 







DQ12 

A 

74^i 

98-ti 

nil 

3-ti 


5-ti 

8-a 


B 

189-ti 

25l^i 

Hi 




9^i 

DQ13 

A 


6-ti 







B 


36-ti 






DQ14 

A 

18^i 


6-ti 






B 

108^,i 


36-ti 



' 
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TABLE III (Contd. ) 


Line 

Section 

Add 

Subtract 

Multiply 

Divide 

SQRT 

Trans . 
Function 

If 

(Branch) 

Initialize 

DQ15 

A 

6^i 








B 

36^i 







DIFEQ loop used to evaluate DQ3, DQ6, DQ9, DQ12 


**7a 

A 

33k^i 

30k-ti 

5k-ti 


9k-ti 

4k-ti 

5k-Li 


A 

Ik-ti 

Ik-ti 






c 

A 

2k<e,i 

2k-tt 






d 

A 

lk-t.i 


Ik-ti 





e 

A 

Ikti 

2k-ti 

Ik^i 





f 

A 


3k^,i 

Ik^i 





g 

A 

Ik-ti 


Ik-ti 





h 

A 

3k^i 

Ok-ti 






i 

A 

3k^i 

Sk'ti 






j 

A 

2k^i 

3kti 


Ik^i 




k 

A 

3k^,i 







1 

A 

2k-ti 

3k-ti 


Ikt 
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TABLE III (Contd, ) 


Line 

Section 

Add 

Subtract 

Multiply 

Divide 

SQRT 

Trans. 

Function 

If 

(Branch) 

Initialize 

m. 

A 



3k^i 





n 

A 

2k-t,i 

3k>a 






o 

A 

5k-ti 

7k-f.i 

Ik-ti 





P 

A 

4k-ti 

5k-ti 

lk>f,i 

lk<ti 




q 

A 

4k-t,i 

18k-ti 

2k^i 





r 

B 

9k-ti 

35k-ti 

Ik-ti 





S 

B 







9kli 

t 

A 







3kli 

u 

A 

6k>ti 

9k-ti 

Ik-ti 





V 

B 

180k-ti 

2 16k-ti 



• 



w 

A 

Ik-t/i 







X 

A 


— 

1 

! 

• 

Ik-ti 


Solution to Kepler's equation used to determine 7a 

KEl 

A 







Ik-ti j 

KE2 

A 

Ik-^ i 

2k-ti 

Ik-ti 


Ik-ti 
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TABLE IlKContd. ) 


Line 

Section 

Add 

Subtract 

Multiply 

Divide 

■I 

Trans. 

Function 

If 

(Branch) 

Initialize 

KE3 

A 





Ijk-ti 



KE4 

A 





Ijkti 



KE5 

A 

Ijk-ti 

Ijk-ti 






KE6 

A 

2jk^i 

2jk-t-i 






KE7 

A 

2jk-ti 

2jkli 






KE8 

A 

2jk^i 

2jk^i 

Ijk-ti 





KE9 

A 







ljk^.i 

KEIO 

A 

Ijk-ti 







KEll 

A 






Ijk^i 
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Section 
of ESF 

Add 

Subtract 

Multiply 

Divide 

SQRT 

Trans. 

Function 

If 

(Branch) 

Initialize 

A 

720i 

1248i 

220i 

40 i 

16i 

lOi 

20i 

B 

356i 

469i 

27i 

li 

9i 

4i 

95i 

C 

46i 

60i 

2i 



li 

I2i 


TABLE IV 


Total Number of Operations Required by One Complete 
Execution of the ESF Algorithm in ”i" Iterations 
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Section 
of CSF 

Add 

Subtract 

Multiply 

Divide 

SQRT 

Trans. 

Function 

If 

(Branch) 

Initialize 

A 

3400i 

4l60i 

860i 

I20i 


200i 

420i 

B 

10, 582i 

12, 062i 

1, I26i 



lOi 

756i 

C 

7i 

6i 




li 

42i 


TABLE V 


Total Number of Operations Required by One Complete 
Execution of the CSF Algorithm in "i" Iterations 









Operation 


Relative Execution Time^ 


Addition 1 

Subtraction 1 

Multiplication 2 

Division 4 

Square Root 3 

Transcendental Function 9 

If (Branch) 1 

Initialize 2 


‘ Approximate "add - times” where one add- time = . 66 milliseconds 


TABLE VI 


Relative Execution Time of Operations on the AGO 
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■ 

A 

B 

C 

Total 

ESF 

4, 410i 

1, 680i 

199i 

6, 289i 

CSF 

19, 800i 

40, 8l2i 

ll3i 

60, 725i 


TABLE VII 

Total Number of "Add Times" Required by One Complete 
Execution of the ESF and CSF Algorithms in "i" Iterations 
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A(sec) 

B(sec) 

C(sec) 

Total 

ESF 

2. 91i 

1. Hi 

. 13i 

4. I5i 

CSF 

13. li 

26. 9i 

. li 

40. li 


TABLE VIII 


Estimate of Total Computation Time for the AGC Required by One 
Complete Execution of the ESF and CSF Algorithms in "i" Iterations 
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18 1 



STEP SIZE ( Degrees ) 

109 


25.5 


1 


10 



ERROR IN SIGNIFICANT FIGURES AFTER 36 MEASUREMENTS 




RUN TIME ( SEC) 


RUNTIME Vs ERROR 



I 23456789 )0 

ERROR (SIG. FIG.) 

Figure 4. 3 Effect of Computer Run Time on Integration Error 
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(m) 


Iiesf " -cspl 



Figure 4. 4 Magnitude of the difference between the estimated position of the ESF and the CSF 
for one Monte Carlo run 
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. 0020 - 

(m/sec) 

. 0015 -- 

.OOlO- 

. 0005 - 

ot 


I V “ V I 3, ” 0 

l-ESF -CSF' -d 

100 n. mi circular earth orbit 


O 

10 measurement intervals 



20 “10 60 

Figure 4. 5 Magnitude of the difference between the estimated velocity of the ESF and the 
CSF for one Monte Carlo run. 


0 



114 


1.5 + 


(m) + 


1.0 + 


.5 + 


0 - 

O 


f-ESF -TRUE' 





a , = 0 
— d 

100 n. mi circular earth orbit 

O 

10 measurement intervals 

O 

itial rms position error = 8. 84xl0""m 

ihitial rms velocity error = 8.65 m/sec 

R 9 

easurement variance = 10 m 





MEASUREMENT 

NUMBER 


Figure 4. 6 The difference between the actual position error of both the ESF and 
the CSF for one Monte Carlo run 
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116 


(m) 

3000- ■ 


2000 + 


1000+ 


0 


□ = rms position error prior to update 

+ = rms position error alter update 


ID 


□ 


□ 


a . = 0 
— d 

100 n. mi. circular earth orbit 

O 

10 measurement intervals 

2 

initial rms position error = 8.84 x 10 m 

initial rms velocity error = 8. 65 m/sec 

6 2 

measurement variance = 10 m 


O' 


iO 


O 


++ 


(D 


O 


+ O 


O 


+ + 

+ Q 

{ Q° 

+ . ° O 

tp'*' a 

+ » (P □ 




+' 


20 


+0 


60 


MEASUREMENT NUMBER 

Figure 4. 8 RMS estimated current position error for one Monte Carlo run 
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0 


8 + 

(m/sec) 


•-0 


8 + 


lO 


0 


rms velocity error prior to update 
rms velocity error after update 

^d = 0 

100 n. mi circular earth orbit 

10 measurement intervals 
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Figure 4. 9 RMS estimated current velocity deviation error for one Monte Carlo run 
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Figure 4. 10 RMS estimated position error at the epoch for one Monte Carlo run 
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Figure 4. 11 RMS estimated velocity error at the epoch for one Monte Carlo run 






Figure 4. 13 Magnitude of the velocity deviation for both the ESF and the CSF for one 
Monte Carlo run 
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Figure 4. 15 Magnitude of the difference between the estimated velocity deviation of the ESF 


and the CSF for one Monte Carlo run 
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* Figure 4. 17 Magnitude of the difference between the estimated velocity of the ESF and 


the CSF for one Monte Carlo run 



126 


I A A 

'-ESF ' i^TRUE 



-CSF -TRUE' 



due to 

100 n mL circular earth orbit 

O 

10 measurement intervals 
initial rms position error = 8, 84x10^ m 
initial rms velocity error = 8.65n/sec 
measurement variance = 10® m^ 


MEASUREMENT NUMBER 


Figure 4.19 Magnitude of the Actual Position Error of the CSF for one Monte Carlo Run. 
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Figure 4. 21 RMS estimated current position error for one Monte Carlo run 
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Figure 4. 23 


Magnitude of the Actual Error in the Velocity Estimate of the CSF 
for One Monte Carlo Run 
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Figure 4. 25 RMS estimated current velocity error for one Monte Carlo run 
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Figure 4. 26 Magnitude of the change in the epoch position from its true initial value for 
one Monte Carlo run 
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Figure 4. 28 Magnitude of the difference between the estimated position of the 
ESF and the CSF for one Monte Carlo run 
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Figure 4. 29 Magnitude of the difference between the estimated velocity of the ESF 
and the CSF for one Monte Carlo run 
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Figure 4. 30 The difference between the actual position error of both the ESF and the CSF 
for one Monte Carlo run 
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Figure 4. 31 The difference between the actual velocity error of both the ESF and 
the CSF for one Monte Carlo run 
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Figure 4. 33 Magnitude of the velocity deviation for both the ESF and the CSF 
for a different Monte Carlo run 
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Figure 4. 34 Magnitude of the difference between the estimated position deviation of the 
ESF and the CSF for a different Monte Carlo run 
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Figure 4. 35 Magnitude of the difference between the estimated velocity deviation of 
the ESF and the CSF for a different Monte Carlo run 
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Figure 4. 36 Magnitude of the difference between the estimated position of the ESF and 
the CSF for a different Monte Carlo run 
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Figure 4. 37 Magnitude of the difference between the estimated velocity of the ESF 
and the CSF for a different Monte Carlo run 
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Figure 4. 38 Magnitude of the actual error in the position estimate of the ESF for 
a different Monte Carlo run 
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Figure 4. 39 Magnitude of the actual error in the position error of the CSF 
for a different Monte Carlo run 
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Figure 4. 41 RMS estimated current position error for a different Monte Carlo run 
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Figure 4. 42 Magnitude of the actual error in the velocity estimate of the ESF for 
a different Monte Carlo run 
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Figure 4. 43 Magnitude of the actual error in the velocity estimate of the CSF 
for a different Monte Carlo run 
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Figure 4. 44 The difference between the actual velocity error of both the ESF 
and the CSF for a different Monte Carlo run 
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Figure 4. 45 RMS estimated current velocity error for a different Monte Carlo run 
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Figure 4. 47 Magnitude of the change in the epoch velocity from its ture initial 
value for a different Monte Carlo run 
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Fig. 4. 49 The difference between the actual position error of both the ESF and the CSF 
one Monte Carlo run 
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Figure 4. 50 RMS estimated current position error for one Monte Carlo run 
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Figure 4. 53 RMS estimated current velocity deviation error for one Monte Carlo run 
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Figure 4. 54 Magnitude of the position 
for one Monte Carlo run 
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Figure 4. 55 Magnitude of the velocity deviation for both the ESP and the CSF for 
one Monte Carlo run 
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Figure 4. 56 Magnitude of the difference between the estimated position deviation 
the ESF and the CSF for one Monte Carlo run 
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Figure 4. 57 Magnitude of the difference between the estimated velocity deviation 
of the ESF and the CSF for one Monte Carlo run 
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Figure 4. 60 Magnitude of the actual error in the position error of the CSF for one 
Monte Carlo run 
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Figure 4. 62 RMS estimated current position error for one Monte Carlo run 
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Figure 4. 63 Magnitude of the difference between the estimated velocity of the ESF 
and the CSF for one Monte Carlo run 
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Figure 4. 64 Magnitude of the actual error in the velocity estimate of the ESF 
for one Monte Carlo run 




Figure 4. 65 Magnitude of the actual error in the velocity estimate of the 
CSF for one Monte Carlo run 
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Figure 4. 66 The difference between the actual velocity error of both, the ESF and the 
CSF for one Monte Carlo run 
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Figure 4. 67 RMS estimated current velocity error for one Monte Carlo run 
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Figure 4. 69 Magnitude of the change in the epoch velocity from its true initial value 
for one Monte Carlo run 


CHAPTER V 


CONCLUSIONS 


From the computer simulation results of Chapter IV it is seen 
that the epoch state filter produces the estimated state of the space- 
craft with nearly the same accuracy as the current state filter, but 
with a considerable savings in computational time. The inaccuracy 
in the ESF, incurred because of the assumption that allows for the 
analytical calculation of the state transition matrix, is small for a 
100 nautical mile circular earth orbit with a disturbing acceleration 
due to as much as 10 J 2 . At times in the orbit, the ESF is even 
shown to be more accurate than the CSF. 

Computer simulation results were for low earth orbits where 
the effect of the disturbing acceleration due to the J 2 term is greatest, 
thus contributing to the largest possible inaccuracy for the ESF. This 
inaccuracy decreases as the spacecraft's altitude increases away 
from the influence of disturbing bodies. It was further demonstrated 
that distizrbing accelerations due to as much as 10 J 2 were tolerable. 
The fact that for a larger circular earth orbit with radius twice the 
equatorial radius of the earth, the errors in the estimate of the 
state of the spacecraft are decreased considerably is consistent 
with this trend. 

The savings in computation time for the ESF on the AGC is a 
substantial improvement over the conventional solution to the 


178 



navigational problem. This savings is due largely to the analytical 
extrapolation of the state' transition matrix as opposed to the numerical 
integration of a matrix differential equation. Also, the integration 
technique for the ESF, which was thought to be more time consuming 
than that of the CSF, proved to be comparable on an IBM 360 Model 
75 computer and faster on the AGC, 

The results of the study indicate that the epoch state filter is 
an economical filter which may be used to estimate the same quantities 
as the Apollo navigation filter. 
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APPENDIX A 


DIFFERENTIAL EQUATION FOR THE 
EXTRAPOLATED COVARIANCE MATRIX 


$ (tj^. y E(t 

k-l* '*''''**k; Vl* 

(A. 1) 

«<‘k- 'k-i' = 

‘k-i> 

(A. 2) 

i^(tk. Vj) = 


(A. 3) 


Differentiating A. 1 with respect to time yields 

+ Vi> 

Substituting Equations A. 2 and A. 3 in A. 4 

+ [«(tk, Vi> E(t^_ j) $’■(1^. j)l F(t^)'T 
Noticing that the terms in brackets are E'(t^) Equation A. 5 reduces to 
E' (tk) = F (y E'(tj^) + E'(y (A. 6) 
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APPENDIX B 


VARIATIONS OF THE 
EPOCH FORMULATION 


Rather than using Q as the independent variable and solving the 
following differential equation 


— = — + [ sin 6 h X r + h (1 - cos 6) r] • a . 

d t r2 h^ 


(B. 1) 


The generalized anomoly, x, may be used instead. The differential 
equation to be integrated is then 


= //M ^ 

U„(x; a) 

V • a _i 

d t r 

_ _d 
M 


where U2(x; a) is the transcendental function 


(B. 2) 


U (x; a) = x"(— - 
” \n! (n 


, 2,2 

(ax ) 


(n+2)! (n+4)! 


-•••) 


(B. 3) 


and 0 / is obtained from the following equation 



(B.4) 


as given in References 1 and 2. Integrating Equation B. 2 for x eliminates 
the need to solve Kepler's equation. This is also the case when B is 
used as the independent variable since Kepler's equation can be solved 
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for either x or 0, The U functions used in solving Kepler's equation 
can be expressed in terms of either x or 0, i. e. , U = U (x; a) or 
U = U(6; a). Still another variational parameter which may be used 
is the variable epoch time, t^, and the equation to be integrated is 

dt - 

— 2 = — (cv+U2r>- a^ (B, 5) 

dt ^ 

where c is defined by 

y^Tc = 3 Ug - xU^ - Ug A/M~(t - t^) (B.6) 

where U2J and Ug are transendental functions of order 2, 4, and 5 
respectively and t is the current time. Use of the differential equation 
for t^ necessitates solving Kepler's equation 


(t - t^) = r^ Uj(x; Of) + Oq U2 (x; a) + Ug (x; a) 


(B. 7) 


for X, In Equation B, 7, is given by 



V n 


V 

— o 


(B.8) 
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APPENDIX C 


DERIVATION OF THE VARIATIONAL 
EQUATIONS FOR - G. -F^, AND F 


The equation for - G in terms of the U functions is given as 

_G = - — U, + -21 U„ (C. 1) 

Differentiating this equation according to the formal rules of differen- 
tiation given in Section 3. 4 results in 

-dG 
d t 


r 


dU. 


da 


d t ./F d t 


^ 2 ^ 


_2.11 
d t 


(C.2) 


where the variational equations for the transcendental functions are 
given by 


dU 


1 _ TT / d ^ _j_ 1 -y d O' 


d t 


= u 


° ^ dt 2 ^ dt ^ 2 




u. 


O' 


dt 


(C.3) 


= U. + — U J - J_ 

dt ^ ^ dt 2 dt ^ 2 dt 


(C.4) 


and that of a by 


d oi 
dt 




I* 3d 


(C. 5) 
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Substituting Equations C. 3 to C. 5 into Equation C. 2 and collecting 
terms yields 


dG 
d t 






^d ^2 


The first term in brackets is F and the second bracketed term is G. 
Making these substitutions in Equation C. 6 results in the following 
variational equation for - G 


dG 

dt 



3d 



/ ^ 1 

V — — 

^ dt 2 


jF+J-r 
dt ^ fi ~ 


3d^ 


Similarly, differentiating the equation for F^ where 


Ff - 

r r 


yields 


'‘“i jw 


dt 


dr 


U, 


1 

rr dt rr dt 

o o 


dU 


Upon substituting Equation C. 3 for 


1 : <z £ da 


and C. 5 for and 

d t d t 


dr 


d t 


= - <7, 


d t 


d t 


-li + J_ u„ ) - J. (r + r^) Up - U, 


a 


d t 


d t 


(C.6) 


"2 

(C. 7) 
(C.8) 

(C.9) 


(C. 10) 
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where 


dg 
d t 



(C. 11) 


Into Equation C. 9, the resulting equation is 



Collecting terms and simplifying Equation C, 12 produces 



(C. 13) 


The second and fourth terms on the right-hand side of Equation C, 13 
cancel and the term in brackets is G^. Making these changes and col- 
lecting related terms results in 
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, JIT / dg + J_ u 

dt r ^ dt 2 ^dt 


N ^1/^2 


r 

o o 


I ■ *d ■" 3 • 5d 


- V • a , + 
— — d 


Uj ./jr 


r r 


!1 • 5d) 


(C. 14) 


where V • has been added to and subtracted from the second term 
in parantheses. Collecting terms within the second set of parentheses 
results in - V • + (-F^ r* + F v* a^) which is equal to (v^ - v) a^. 

Thus, the variational equation for -F^ is 


dF 


i = flA + _l 


d t 


I- 


U, 


dt 


d a 
dt 


G, - ^ (v^ 


v) 


3d^t 


(C. 15) 


Finally, the equation for F is 

F = 1 - (C. 16) 

r 

o 

Differentiating this equation as was done previously for -G and F^ 
yields 


Substituting for 


dF 
d t 


dU„ dr 

2 O 


dU 
d t 


2 ^ "2 "'■o 


d t 


(C. 17) 


dt dt 

C. 5, and C. 11 respectively results in 


dcy d (7 

— , and with equations C. 4, C. 10, 

dt dt 
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dF 
d t 


1 






. Vr . UjU 2 r-a^ 

+ 2 ^ ^d-^— -^d-— 2 


Mr 


Mr 


JT 


(C. 18) 


Cancelling the second and fourth terms on the right-hand side of 
Equation C. 18, replacing the term in brackets by ^n^ G, and collect- 
ing related terms yields 


dF ^ 
dt 


±. (AL + ± xj l£L ) G - 

r„ ^ d t 2 d t ^ 
o 


rU2 ,V JW 
I r • a 

nr r r ~ 

^ o o 


^2 

V • ^d V* a^ - V- a^) (C. 19) 


where v • a^ has been added to and subtracted from the term in the 
second set of parentheses. But this term is just — (v - v ) a . so 

jU -o - -d 

that the variational equation for F is 


dF 

dt 


^ (-ii + — U, — ) G + i (v - v) • 
r 2 Wt 2 3 dt ^ /J - 

o ^ 




- JL (v - v) . a . F 
^ -o - -d 


(C. 20) 
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EPOCH STATE FILTER 


BEGIN 

? 
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INDEX I,, I , 

DIM6NSI 0M( E , 6X6) ,( B » 6X 3 ) , ( PH I , 6X6 ) , ( W , 6 ) , ( RD , 6 ) , ( OXF, 6 ) , 

{ DO,BO) ,( nxOF.6) ,( PHI 1,6X6) ,(PHI2,6X6),(PHI3,6X6),(DVTM,6) 

nn in ? for i = o( i )?9 

DO = RNDMN( 1000) ' 

I 

E = 0 

no TO 3 FOR I = 0( 7 ) 14 • 

E = 26OTO0 
I 

DO TO 4 FOR I = 21( 7 )35 

E •= 

I 

■f 

R = 0 

DO TO 6 FOR I = 0( 4)8 

R = 1 
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1P5P1 = 1 .6 PI 
3PI = 3 PI 
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alpha =1-010 

15 

MU = .398601? 10 
RMU = SORT! MU ) 

RR = 6563365 

R = RR . ■ 

VC = SORT! MU / RR ) 
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PDOSC = ? PI 50RT( RR / MU) 

DT = PDOSC / 7? 


RE = 6378165 
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J2 = .noiop23 
J = 0 
T = n 


THETA 

= n 



thetai = 

0 


IZ = 

( 0. 

(S 

3 3 

RC) = 

( RR , 

0, 

0) 

ROI = 

RCI 



IRO = 

ONT T( 

RO I 

RO = 

RR 



R . =• R 0 



IR = 

1 RO 



VO = 

( 0, 

vn. 

(33 

vni = 

vn 



OROE 

= ( 0 

, 0 

. 0) 

ORE = 

DROE 


DVOE 

= DRUE 


DO TO 

1 6 

FOR 

T = 


ITH = IZ=MR 
' ITHQ .= IZ’^IRD 

DELTH = ( ITH.DRE/R - ITHO.DROE/RO ) 
THETA = THETA + DELTH 
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? 7. 7 - - 

AD = ( -Mll/R ) 1.5 ^7 ( RF/R ) ((1- 5 CPHI ) IR + ? CPHI 17) 

CTHETAI = r,DS( THFTAl) 

STHFTM = FUKThETAD 

ROI = ARVAI ( RPT ) 

vni = arvai ( vni ) 

SIGMAHl = RDl.Vni / RI^K 

2 

ALFA! = ( ' / ROI ) - ( Vni / I'll! ) 

2 2 

]P = ? RUT - At. EAI ROI - SIGH AO I 
HI = SORT! P'll IP ) 

RI = IP/ (1+ ( IP/ROI -1) CTHETAI - HI SlGMAOl ETHETAI/ ROI RI'^U) 
FI = I - PI (1 - CTHETAI )/ IP 

ETI = ROM EKIMAOI (1 - CTHETAI)/ RCU IP - MU STHETAl / ROI HI 

f;i = RI RMI RThFTAj / HI 

GTI = 1 - ROI ( 1 - CTHETAI ) / IP 

RI = >1 RMI + Gl vni 

IRI = MMIT( RI ) 

CPHI 1 =■ IRI . 17 

VI = ETI KOI +GTI vni 

HI = RI=;<VI 

- 2 2 2 - 
AnI = ((-MU/RI ) 1.5 J2 (RE/Rl) ((1-5 CPHII ) IRI + 

2 CPHI I TZ) ) 



E 

M 

E 

M 

M 

M 

M 

M 

E 

M 

e 

M 

E 

H 

E 

M 


E 

M 

M 

M 

M 

M 

M 

E 

M 

E 

M 

E 

M 

E 

M 

E 
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RO = Rl) + DROE 

VO = VO + DVOE 

DO TO 13 EOR L = 1( 1)2 

no TO 1 3 FOR K = 0( 1 )3 ■ 

CTHET/i = r,OS( THETA) 

STHFTA = '^INK THETA) 

RO = APVAi ( RO ) 

VO = AHVAi, ( vn ) . • 

SIOMAO = KO.VG /,RMU 

2 

ALFA = ( ?■ / RO) - ( vn / MU) 

2 2 

p = ? RO - ALFA RO - SIGMAO 
. H = SORT) MU P ) . 

R = p / (1 + (p/RO -1) CTHETA - H SIGMAO STHETA / RO RMU) 

F = 1 - R ( 1 - CTHETA ) / P 

FT = RMU '^IGMAO (1 - CTHETA)/ RO P -.MU STHETA / RO H 

G = R RO "CTHETA / H 

GT = 1 - RO ( 1 - CTHETA) / P 

R = F RO + G VO 

IR = UNIT) R ) 

CPHI = IR.IZ 

V = FT RO + GT VO 

H = R>:'V 


5 



2 - - ?■ 

OTHETA/DT = H / R + ( STHETA H>i‘R + H (1 - CTHETA) R ).AD / H 

DRO/DT = (R Rn/ (MU P)) (1- CTHETAJ ( ( RO- R) V + V R ) AD - G AD 

DVQ/DT = ( R / MU) ( VO - V) ( VO - V ) AD + F AD 

? - - - 
DTHETAI/OT =(HI/R1 + (STHFTAl HI^fRI + HI (1 - CTHFTAI) RI). 

? 

ADI / HI ) 

DROI /rn =((RI ROW (MU IP)) (1 - CTHFTAI) ((ROI - RI) VI + VI 
R I ) ADI - GI Ani ) 

OVni/DT = (RI / MU) (VOI - VI) (vni - vn ADI + FI ADI 
OIFEO T,DT .OTHFTA/DT,nRO/nT,DVO/DT ,DTHETAI /DT.DROI /DT ,DVOI /OT 


ROC MG = AKVAI ( 

RO 

- ROI ) 


_ 

_ 

VnCMG = AHWAI ( 

VO 

- VOI 

RO = AHVAi ( RO 

) 


IRO = IIMIK RO 

) 



IF THFTA <= IP^PI, GO TO 14 
IF THFTA GORFO 3PI» GO TO 12 

CALL KWK.r.OMICSt 4, 0, RO, Vl.U lPSPl, MU 

RESUME TF). RI, VI, THETAl, SLR, COGA, SMA, PERIOD, RI, SOLMTAG 

CALL KWK.i'OMlCS, 4, 0, RI, VI, ( THFTA - IP5PI ), MU 

RESUME TF'.’, R, V, THFTA2,’ SLR, COGA, SMA, PERIOD, R, SOLMTAG 
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CALL REP.rnMICS, 3, 0, RO, VO, TFT, MU, 0 

RESUME TFT.Rl ,V1,THFTAI,SLR ,COG A , SM A , PF R I OD , R 1 , SOL NT A C , ANOMDL Y , 
PHI1 

CALL PEP.'OMICS, 3, 0, Rl, VI, TF?, MU, 0 

RESUME TE ,R , V , THET A? , SL R , COG A , SM A , P ER 1 OD , R ,SDL NT AG , AMOMOL Y , 
PHIP 

PHI = PHI,- PHTl. 

GO TO IS 

IP CALL KWK.r.uwir.s, a, o, ro, vn, ipspt, mu 

RESUME TFI, RT, VI, THETAl, SLR, COGA, SMA, PERIOD, Rl, SOI MTA(; 
CALL KWK.r.fiMICS, 4, 0, Rl, VI, 1P3PI, MU 

RESUMF TEP, RP, VP, THFTAP, SIR, COGA, SMA, PEftlOD, RP, SOL MTAG 
CALL KWK .COMICS, 4, 0, RP, V?, (THETA - 3PI), MU 
RESl.JMF TE3, r, V, THFTA3, SLR, COGA, SMA, PFRIOU, R, sniMTAG 
CALL RFp.rOMir.S, 3, 0, RO, VO, TFT, MU, 0 

RESUMF TFI,R1.V1 , THET AT, SLR, COG A, S(M A , P F R I 00 , R 1 , SOL MT AG , AMOMOL Y , 
* 

PHIl 

CALL REP. CONICS, 3, 0, Rli VI, TFP, MU, 0 . 

RESUME TFP,R?.VP,THETAP,SLR,COGA,SMA,PFRIOD,R?,SOLNTAG,AMnMnLY, 
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CALL PFP.CnNICS, 3, 0, R2, V2, TF3, MU, 0 

RESUME TF3 ,R , V , THETA3 , SL R , COG A , SM A , PER I OD , R , SOL NTAG , AMOMOL Y , 

PHI3 

>!« 

PHI = PHI3 PHI? PHIl 
GO TO 15 

14 CALL PWK .(.nMics, 4, n, Ro, vn. theta, mu 

RESUME TF. R. V, THETA?, SLR, CPGA, SMA, PERIOD, R, SULNTAG 
CALL REP.iOMlCS, 3, 0, RO, VO, TF, MU, 0 

RESUME TF ,R , V .THETA ?, SLR, COGA, SMA.PFR IDO, R , SOL NT AG , AMOMOL Y ,PHI 
T -C 

15 BO = PHI B 

.1 

IR = ,UMIT( R ) 

AO = BU.( I- RO ) + AL PH A 
W = E B() ./ All 

>|s :|£ _ _ 

g = E - W PO E 


DVTM 

= Rn 

- ROI 

0 

0 

0 

DVTN 

= RU 

- RDI 

I 

1 

1 

DVTM 

= RO 

- ROI 

? 

? 

? 

DVTN 

= VI) 

- VO I 

3 

0 

0 



|V| 

DVTM. 

= VO - vni 

S 

4 

1 1 

M 

DVTN 

= VII - vni 

S 

5 

? ? 

E 

M 

nxoF = 

W ( DO 

S 


1-1 

M 

DROE 

= DXOF 

S 

0 

( ) 

M 

DRriF 

= DXOE 

s 

1. 

1 

M 

DROE 

= . nxOE 

S 

? 



Dvoe 

= nxriE 

s 

0 


M 

DVOE 

= DXOE 

s 

1 

4 

M 

DVOE 

= DXOE 

•S 

? 


E 

_ 

— 

M 

DXE =, 

PHI DXOF 

M 

ORE = 

DXF 

s 

0 

0 

M 

DRE - 

DXE 

S 

1 , 

1 

M 

ORE = 

DXF 

s 

2 

,? 

M 

DVE = 

DXF 

s 

0 

3 

M 

DVE = 

DXE 

S 

1 

4 

M 

DVE = 

DXF 

S 

7 . 

5 

M 

J = J 

+ 1 


M 16. I F J > ? , J = 0 


PO.nVTN ) 
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START- AT RFGIM 



CONVENTIONAL STATE FILTER 
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BEGIN 

? 


4 


5 


6 


INDEX T_. .1, 1 

D I MENS I r)M( F,6X6),(B,6X3),(W,ft),(E,6X6),(in,6X6),(PHI,6XB), 
( no. HO) ,( nXF.f>) ,( DVTM»6) 

nn Td ? FOR I = o( 1 ) 79 

DO = RNONM( 1000) ■ 

I 

E = 0 

DO TO 3 FiiR I =0(7)14 

6 = Pf-,0100 

I 

nn Tfl 4 FOR I = ?1(7)35 

F = ?s 
I 

in = o 

00 TO S FOR T = 0(7)35 

in =1 
I 

■f- 

R = 0 

DO TO 6 fur I = 0(4)H . 

R = I 
I 

F = 0 

F = 1 

3 

F = 1 

10 - 

F = 1 

17 

J2 = .0010823 

6 

ALPHA = 1.0 10 



RE 


637ft I 65 


15 

MU = .3936017 10 
RR = 6563365 

3 

PDGSC = 7 PI SORTlRR 7 MU) 
DT = pnnsr. / 360 
VC = SORT! MU / RR ) 

T1 = 0 
vi = n 


IZ = 

( 0, 

0, 

1 ) 


MU = 

( 0, 

0, 

0) 


DEL TA 

= 

NU 



DEL AD 

= 

Mil 



NUAD 

= NU 



RU = 

( RR 

t 0. 

0) 


vn = 

( 0, 

VC , 

0) 


DRE = 

( 0 

. 0. 

0) 


DVE = 

ORE 



Dn TO 

1 ft 

FOR 

I = 

1( 1 )80 

PHI = 

n 




DO TO 

1 1 

FOR 

Z = 

0( 7 )35 


PHI - 1 
Z 



F 


r 

M 

DELTA 

= liFi TA + DRE . 


E 

M 

NU = MU + nvF 

0(1 T(J 1 6 FDR 1 = 1(1)10 


|V( 

nn TO 

1 4 for K = 0{ 1 ) 3 





'' 

M 

ANGL F 

= 2 Pi Ti / pnnsc 


E 

M 

Rose 1 

= ( RR r,ns( ANGLE), 

RR SIN( ANGLE ) , D ) 

F 

|V| 

R = ROSCl + DELTA 


E 

M 

IR = UIMTT( R ) 


F 




H 

R = ARVAl. ( R ) 


6 




Fi 

CPHI = 

IK. 17 


E 


- _ — 

2 

h 

0 = ( ( 

0F( TA - 2 R ). DELTA ) / R 

E 


2 

1.5 

jvi 

FO = 0 

( 3 + 30 + 0)/ 

( 1 + ( 1 + 0 ) ) 

E 

- 

2 

2 2 - 

M 

AD = ( 

-MO/R, ) 1.5 J2 ( 

RE/R ) ((1-5 r.PHI ) IR + ? CPHl I/) 

F 

=;= 

s - - 

2 

1^1 

t; = ( 

HLI / R ) { 3 R R 

- R 1 DMATR IX ) 

M 

F 



s 

1 R 

0 


M 

F 

G 


s 

19 
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M 

F 
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S 

?n 
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H 

■ F = 

G 


S 

24 

3 
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G 
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4 


M ' 

F = 

G 


,S 

26 

5 








M 


F = G 



S 


30 f 



M 


F = 0 



s 


31 7 



M 


F = 0 



S 


3? 8 



E 
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RT = RDSCl + OELAn 



E 


_ _ 



M 


IRT = IIMITI RT ) 



F 


- 



C'l 


RT = AhVAI. ( RT ) 



F 


- _ 



M 


TCPHI = IRT.IZ 



E 

M 


- - - 2 
OT = ((OElAD - ? RT),nFLAO) / RT 


e 


2 


1.5 

M 


FCIT = OT ( 3 + 3 OT + OT ) 

/ ( 1 + ( 1 + 

OT ) 

E 


? 

? 

? 

M ' 


ADT = (-MII/RT ) 1.5 J? (RF/RT) ((1-5 

TCPHI 

E 


_ _ 



M 


DOELTA/OT = Mil 



E 

M 


3 

niMII/DT = ( -Ml) / RR ) ( 

FO R + ORLTA 

) + AO 

F 


_ _ 



M 


DOELan/nT = imuao 



E 

M 


- 3 

DNIIAn/nT = ( -MU/RR ) (FOT 

RT + OFI. AO) 

+ AOT 

E 

M 


s|c >:c 

nPH)/riT = F pht 



E 

M 

] 6 

niFEO T1 , OT. DDELTA/OT, 

DMi l/DT, on EL 

AD/DT, 

E 


. _ 



M 


VOSC = ( -VC STM( ancle ) 

, VC cns( ANGLE ), 

e 

M 


VT = VUSC + M|)AD 




. V = v(j-';c + Mil 


) 

) IkT + ? 


OMDAn/nT , 
0 ) 


TCPHl 17) 


l.'PH] /[)■{ 
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* sit T 

E = PHI F PHI 

-r. * 

A = fi . ( E B ) + ALPHA 

J vt 

— — C 

W = E H ' / A 

J 

=1= - -C ^ * - t-C T - - 

E=(TD-WR )E(IO-WR ) + ALPHA W W 

,1 J 


dvtm 

= DPITA 

- DELAD 

0 


0 0 

DVTN 

= DELTA 

- DELAD 

1 


1 1 

ovtn 

= DELTA 

- DELAD 

? 


2 2 

DVTn 

= Ml 1 

NUAD 

3 

0 

0 

ovtn 

= NU 

NUAD 

A 

1 

1 

DVTn 

= Nil 

NUAD 

5 

2 

2 

_ 

_ 

-C - 

Dxe = 

W ( DO 

- R .DVTN 


I 

-1, . J 

DRE 

= DXE 


0 

0 


DRE 

= DXE 


1 

I 


DRE 

= DXF 


2 

2 


DVE 

= DXE 


0 

3 


DVE 

= DXF 


1 

4 


DVE 

= DXE 


2 

5 


J = J 

+ 1 




18 I F J, > ? . >1 = 0 

CALL viWH.r;RAPH( SII8R ) ,1,2.0, 1 ,5,80,1,0,21,2860,0 


START AT BFGTM 



